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We present the exact solution of the one-dimensional extended Hubbard model in the atomic limit 
within the Green's function and equation of motion formalism. We provide a comprehensive and 
systematic analysis of the model by considering all the relevant response and correlation functions 
as well as thermodynamic quantities in the whole parameter space. At zero temperature we identify 
four phases in the plane ([/, n) (U is the on-site potential and n is the filling) and relative phase 
transitions as well as different types of charge ordering. These features are endorsed by investigating 
at T = the chemical potential and pertinent local correlators, the particle and double occupancy 
correlation functions, the entropy, and by studying the behavior in the limit T — » of the charge 
and spin susceptibilities. A detailed study of the thermodynamic quantities is also presented at 
finite temperature. This study evidences that a finite-range order persists for a wide range of the 
temperature, as shown by the behavior of the correlation functions and by the two-peak structure 
exhibited by the charge susceptibility and by the entropy. Moreover, the equation of motion for- 
malism, together with the use of composite operators, allows us to exactly determine the set of 
elementary excitations. As a result, the density of states can be determined exactly and a detailed 
analysis of the specific heat allows for identifying the excitations and for ascribing its two-peak 
structure to a redistribution of the charge density. 



I. INTRODUCTION 



It goes without saying that the Hubbard model [l| is 
of seminal importance in the study of modern condensed 
matter theory. Originally introduced to describe the cor- 
relations of electrons in the narrow d band of transition 
metals, the Hubbard model has been subsequently used 
to investigate a large variety of situations experimentally 
observed. Upon varying the model parameters, it is be- 
lieved that the Hubbard model can describe many prop- 
erties of strongly correlated electronic systems. However, 
its applicability to real materials is rather limited and 
some extensions are needed. Among various generaliza- 
tions, one of the most studied is the so-called extended 
Hubbard model (EHM), where a nearest- neighbor inter- 
action V is added to the original Hubbard Hamiltonian, 
which includes only an on-site interaction U . If one as- 
sumes that the interaction parameters arise solely from 
Coulomb interactions, then both U and V are repulsive. 
However, these parameters may represent effective inter- 
actions. Thus their ranges can be broader, allowing for 
both positive and negative values. The electronic proper- 
ties of the model are substantially modified by including 
nonlocal interactions allowing one to study general fea- 
tures of doped systems like copper-oxide materials JJ, 
high-temperature superconductors Q, manganese com- 
pounds [J], and to describe charge ordering (CO), exper- 
imentally observed, in a variety of systems. Indeed, sev- 
eral studies of the extended Hubbard model havepointed 
to the presence of states with phase separation [5| and of 



charge-ordered phases [6|,|7|. Furthermore, the introduc- 
tion of an intersite interaction can mime longer-ranged 
Coulomb interactions needed to describe effects observed 
in conducting polymers [8]]. 

Although the Hubbard model is oversimplified and in 
spite of very intensive study, both analytical and numer- 
ical, exact solutions exist only in one or infinite dimen- 
sions. In ID the thermodynamics of the Hubbard model 
has been widely investigated by means of different exact 
approaches for both nearest-neighbor hopping [^, UOl and 
long-range hopping 11]. For higher dimensions {D > 1), 
the study is far from being completed and only few rig- 
orous results are known. In particular, relevant results 
were obtained by exact diagonalization on small clusters 
[l2| and by numerical investigations [l3| • For the case of 
infinite dimension, an exact solution has been obtained 
by dynamical field theory |14| . 

The situation is more severe for the extended version 
of the model: even in ID there is not a complete exact so- 
lution. A particular extended version has been examined 
in Rcf. [l5| , where a kinetic term (Y) for paired carriers 
is added to the usual electron hopping (t). The model 
has been solved only in the particular case of < = —Y. 
Actually, it is known that the EHM is non-integrable 
for general values of the model parameters [10|. These 
difficulties have led to the investigation of the so-called 
narrow band limit, where the hopping t can be neglected 
with respect to the interaction couplings U and V. This 
simplified version, also known as the atomic limit of the 
extended Hubbard model (AL-EHM), has been largely 
studied for several years. 



There is a large and diffused interest in the study of 
this relatively simple model, motivated by several rea- 
sons; (i) the phase diagram has a very rich structure 
with a multicritical behavior; (ii) one observes several 
kinds of long-range orders pertinent to commensurate 
sublattice charge-density alternation (long-range order 
has been observed in many narrow-band materials, such 
as BaBi^Pbi-xOs, Cs2SbCl(i, Fe^Oi, CsAuCls, and 
several tetracyanoquinodimethanide salts); (iii) the re- 
sults for the atomic model may give some information 
beyond the zero-bandwidth limit [several approximate 
studies, as well as some exact results for the ID case 
[l7| . hint at the possibility that charge ordering persists 
in the corresponding narrow band model] ; (iv) the knowl- 
edge of the atomic limit can be used as a starting point 
for a perturbation expansion in powers of the hopping t. 

Without pretending to be exhaustive, a review of the 
different approaches available in the literature could be 
presented as follows. The phase diagram of the AL-EHM 
at half filling {n = 1) has been investigated in several pa- 
pers by means of a modified Hartree-Fock theory [18| . 
beyond Hartree-Fock approximation employing a decou- 
pling procedure [19| , variational approach [20| and Bethe- 
Peierls approximation [21|. The molecular field approx- 
imation [22i] has been used to obtain the phase diagram 
in the (n, U/zV) plane and in the (/x, U/zV) plane {z be- 
ing the coordination number and /i the chemical poten- 
tial) . Another approach makes use of the transfer-matrix 
method that, although practically employable only for 
one-dimensional systems, can provide exact results. In 
Refs. [23, [241 the specific heat and the static magnetic 
susceptibility have been studied at half filling. The case 
of quarter filling, in the limit U — > 00, has been consid- 
ered in Ref. 2J] , where exact expressions for the charge- 
charge correlation function have been obtained for n — 1 
and n = 0.5 {U -^ cx)). For n = 1 and U > 2V, the 
ground state is a Mott insulator with one electron per 
atom; when U < 2V the ground state is a CO state with 
one sublattice doubly occupied and the other empty. The 
phase transition (PT) from the CO state to the Mott 
state can be either first or second order, and the phase 
diagram exhibits a tricritical point. 

A method of constructing the ground state phase dia- 
gram has been formulated in Refs. [2^ by exploiting the 
reflection positivity property with respect to reflections in 
lattice planes exhibited by the AL-EHM. Estimates of the 
probability of some events at nonzero temperatures were 
also provided. In Refs. [2g| the AL-EHM has been stud- 
ied in 2D addressing the problem of determining the zero 
temperature phase diagram in the {fJ^/V, U/V) plane. 
The ground state exhibits several kinds of charge long- 
range orders, which persist also at sufficiently low tem- 
peratures, provided one stays sufficiently far away from 
the boundary between the different regions. This study 
reports the absence of PTs among the different regimes 
at finite temperatures, in contrast to the results of the 
mean- field approximation [2^. It is also suggested that 
the staggered charge order persists in the corresponding 



narrow band model, although there is no rigorous proof 
since the reflection positivity property fails to be true for 
non-zero t. 

Another category of studies is based on the use of the 
Pirogov and Sinai methods [27| . This approach, origi- 
nally conceived to investigate classical spin systems at 
low temperatures, has been extended to quantum sys- 
tems: spin, fermionic and bosonic models. When ap- 
plied to the AL-EHM, these methods allow for a detailed 
study of the low temperature phase diagram. Moreover, 
they allow for an extension to non-zero t, treating the 
narrow band Hubbard model as a quantum perturba- 
tion of the t = model. In Ref. [2^| the phase dia- 
gram of the AL-EHM has been derived. At zero tem- 
perature, the phase diagram in the {U, /x) plane shows 
six phases: three homogeneous (characterized by parti- 
cle densities n — 0,1,2 and a staggered charge density 
A = 0) and three staggered (characterized by particle 
densities n — 1/2,1,3/2 and staggered charge density 
A = 1/2,1,1/2). The staggered order parameter be- 
ing defined as A = limA^oo J2i £(*)'^(*)/^: where e{i) 
assumes two values, 1 or -1, depending on which sub- 
lattice of the bipartite lattice A the site i belongs to. 
For any sufficiently low non-zero temperature, the three 
staggered phases merge into one, without any PT. These 
results have been confirmed in Ref. [2^, where also the 
case of attractive intersite potential {V < 0) has been 
considered. The narrow-band extended Hubbard model 
in 13 > 2 has been studied in Ref. [301 , where it has been 
shown that the staggered charge order persists for suffi- 
ciently small values of the hopping parameter t. A com- 
prehensive study of the phase diagram of the ID EHM 
was carried out in Ref. ^31] using different methods: per- 
turbation analysis in the strong and weak coupling limits, 
quantum Monte Carlo in the intermediate coupling re- 
gion and exact diagonalization calculations. The authors 
individuated various phases, including charge and spin 
density waves, superconductivity, and phase separation. 
G-ology investigations [33, |33|, renormalization group 
[3J] , and bosonization [35| studies have provided analytic 
insights, particularly in the weak coupling regime. 

Finally, a large activity by numerical simulation must 
be mentioned. The EHM has also been studied by means 
of quantum Monte Carlo (QMC) simulations [Sy, [33], 
Lanczos technique |38|, exact diagonalization calcula- 
tions [39|, [40| , density-matrix renormalization group [4l| . 
The phase diagram and some thermodynamical prop- 
erties at half filling of the EHM have been studied by 
means of the density-matrix renormalization group ap- 
plied to transfer matrices |42] . Recently, some interesting 
results in the study of the AL-EHM have been obtained 
by means of Monte Carlo (MC) simulation, for dimen- 
sions greater than 1. In Ref. ,43], by combining MC and 
mean-field approximation, the phase diagram at fixed /x 
has been derived both at T = and finite temperatures, 
and the critical behavior in the vicinity of the tricriti- 
cal points has been analyzed. In Ref. [4j] MC simula- 



tion has been used to study the AL-EHM on a square 



lattice {L — 10-40); a detailed study of the model at 
finite temperature has been presented for different val- 
ues of the filling and/or chemical potential, evidencing 
the presence of different charge orderings and tricritical 
point lines, and the possibility of phase separation and 
formation of stripes. The variational cluster approach 
has also been successfully applied for the EHM in one 
and two dimensions, showing that QMC and density- 
matrix renormalization-group results can be reproduced 
with very good accuracy [4^ . Most of these analyses have 
focused on the half filling case, paying attention to the 
boundary line separating spin-density wave from charge- 
density wave, as well as to the existence and location 
of a tricritical point. In the last few years, there have 
been discussions galore on the existence of an interme- 
diate bond-order wave phase [46| , whose presence seems 
widely confirmed by now [47| . The case of quarter filling 
has been examined in several works, providing precise es- 
timates of the charge gap [4^ , hinting at the possibility 
of superconductivity in an unexpected region of the (C/, 
V) plane [49!| and the existence of a CO transition p, [40| ■ 

As it emerges from the above review, an exact, rigorous 
and detailed study of the EHM and, in particular, of the 
one-dimensional case, for different values of the parame- 
ters n, T, U, and V, is still missing in the literature. Most 
of the existing studies have been restricted to the case of 
half filling. The analysis given in Refs. |26l . |28| . l29l . |30 | 
allows for a complete description of the phase diagram 
at zero temperature, but in the (/i, U) plane. More- 
over, in several works, many properties of the system 
are presented as a function of the chemical potential /i. 
Although from a theoretical point of view the chemical 
potential and the particle density n are conjugated vari- 
ables, experimentally it is possible to tune the density 
in a controlled way, by varying the doping. The chem- 
ical potential is an internal parameter, dynamically de- 
termined by the system itself according the values of the 
external parameters. It is very difficult to induce trans- 
formations at fixed chemical potential. Moreover, by fix- 
ing the particle density or the chemical potential one can 
obtain different results, if the canonical and grand canon- 
ical ensembles are inequivalent due to a divergence of the 
compressibility. For instance, if the chemical potential 
is considered as an external parameter, one finds that 
a solution exists only for discrete values of the particle 
density (n — 1/2, 1, 3/2,2). Thus, all the studies where 
/J, is considered as an independent variable are doubtless 
interesting from a mathematical point of view, but are 
not directly related to experimental situations. It is just 
more physical to consider the particle density number 
n = Ne/N - where N^ is the total number of electrons 
and N is the number of sites - as an external parame- 
ter, determined by the particular material under study. 
In real materials n can vary in a large range: in cuprate 
materials and in conducting polymers, n = 1 (half filling) 
in the insulating state, but largely varies under doping 
[8|; in charge transfer salts n varies by varying the pres- 
sure in the neighborhood of n = 2/3 i50l|; Bechgaard salts 



have n — 1/2 (quarter filling) [5l| . 

It is worthwhile to recall that there is a correspondence 
between the atomic limit of the extended Hubbard model 
and the spin-1 Ising model defined on the same lattice. A 
transformation from the fermionic to the spin Hamilto- 
nian can be performed by defining the pseudospin vari- 
able S{i) = n(i) — 1. Under this transformation, one 
ends up with a Hamiltonian pertinent to a spin-1 Ising 
model with nearest-neighbor exchange J = —V, in the 
presence of an effective crystal field A = U/2 + ksT log 2 
and of an external magnetic field h = ji — 2V — U/2. 
The temperature dependence of the crystal field is due 
to the spin degeneracy induced by the mapping between 
S and n. The literature on the spin-1 Ising model is 
vast and we refer the interested reader to the bibliog- 
raphy given in Ref. |52l |. However, here we would like 
to comment on the equivalence between the two models. 
There is a dramatic difference from a physical point of 
view because, for a spin system, it is possible to tune the 
magnitude of the applied magnetic field h and vary the 
temperature T. The system responds to such variations 
by rearranging the spins in a pertinent configuration, de- 
scribed, in average, by the magnetization m = {S). For 
the fermionic model, the external thermodynamic param- 
eters are n and T: the system responds to variations of 
these parameters by adjusting the chemical potential ^. 

In this work we present a study of the ID AL-EHM 
by using a new method based on the equations of motion 
and Green's function formalism. Recently it has been 
shown [53] that a system built up of q species of local- 
ized Fermi particles subject to finite-range interactions, 
is exactly solvable in any dimension. By exactly solvable 
we mean that it is always possible to find a complete set 
of eigenenergies and eigenoperators of the Hamiltonian 
closing the hierarchy of the equations of motion. Thus, 
one can obtain exact expressions for the relevant Green's 
functions and correlation functions depending on a finite 
set of parameters, whose knowledge is required in order 
to have a complete solution of the system. In the case of 
a one-dimensional lattice, for (7 = 1, 2, 3 152, Hj, [53| , and 
in the case of a Bethe lattice for q = 1 [5g| , it has been 
shown how it is possible to exactly fix such parameters by 
means of algebra constraints. The analytical exact solu- 
tion of the ID AL-EHM was reported in Ref. [52| , where 
preliminary results were given for vanishing on-site po- 
tential. Here we present a systematic and detailed study 
of the model. We study the properties of the system 
as a function of the external parameters n, T/V, U/V 
(throughout the paper we set T^ > as the unit of energy) 
allowing for the on-site interaction U to be both repulsive 
and attractive. Owing to the particle-hole symmetry, it 
is sufficient to explore the interval [0, 1] for the param- 
eter n. The chemical potential is then self-consistently 
determined as a function of the external parameters. We 
also address the problem of determining the zero tem- 
perature phase diagram in the {U/V, n) plane. Relevant 
local quantities, such as the double occupancy, the next- 
nearest correlation functions, the density of states, and 



the correlation functions for the charge and double occu- 
pancy are systematically computed both at T = and as 
functions of the temperature. 

For specific values of the particle density, we find reg- 
ular spatial distributions of the electrons determined by 
the two competing terms in the Hamiltonian, i.e., the 
on-site and the intersite interactions. For particle densi- 
ties less than quarter filling we do not observe CO states 
by varying the on-site potential U . At U — Q there is a 
CO transition at quarter filling from a spatial disordered 
configuration where the electrons are coupled in pairs 
to a checkerboard distribution of singly occupied sites 
for C/ > 0. For particle densities between quarter and 
half filling, we find three possible spatial arrangements 
of the electrons by varying U . For attractive on-site po- 
tential the electrons are still disorderly coupled in pairs; 
at [/ = one observes a transition to a peculiar ordered 
state with alternating empty and occupied sites. The 
peculiarity consists in the fact that, upon fixing the par- 
ticle density, some of the singly occupied sites become 
doubly occupied - with the constraint that the double 
occupancy D is equal to (n — 1/2) - and the possible dis- 
tributions along the chain are energetically equivalent. 
Further increasing the on-site potential, a second transi- 
tion to a disordered distribution of only singly occupied 
sites is observed at t/ = 2V . At half filling, we find the 
well-known transition at C/ = 2V from a state where ev- 
ery second site is doubly occupied to a state with one 
electron per site [18] . We also find the interesting result 
that the correlation functions arc controlled by two cor- 
relation lengths S,n and ^c, though the latter is relevant 
only at T 7^ 0. 

On the other hand, if one fixes the on-site potential and 
allows the particle density to vary, one observes PTs only 
at quarter and half filling, signalled by the divergence of 
the correlation length ^„. For fixed attractive on-site 
potential, the correlation length ^„ diverges only at half 
filling and thus, a long-range order is established. For 
values of on-site potential < C/ < 2V the correlation 
length ^„ diverges for all values of the filling in between 
quarter and half filling, but only at these extremal values 
does one observe a CO state. For large repulsive on-site 
potential, ^„ diverges only at quarter filling, where one 
finds a regular alternating distribution of singly occupied 
sites. The further increase of the particle density leads to 
an increasing single occupation of the empty sites until 
the homogeneous distribution is reached at half filling. 

At finite but low temperatures, one expects that the 
scenario previously depicted still holds. This wisdom is 
strengthened by the behavior of thermodynamic quanti- 
ties such as the compressibility, the charge and spin sus- 
ceptibilities, the internal energy, the specific heat, and 
the entropy. The computation of the entropy requires 
some attention since the ground state is macroscopically 
degenerate (unless for U < 2V and n = 1) and we find, 
thus, a finite zero-temperature entropy. 

Before ending this section, we would like to stress the 
motivations of our work. Firstly, we use a new method. 



The transfer matrix method is the most common ap- 
proach used in the literature to obtain analytical and 
exact solutions of the AL-EHM and, thus, allowing for 
performing explicit calculations of the relevant correla- 
tion functions and thermodynamic properties. However, 
this method is basically restricted to ID systems. There 
are great difficulties in extending the method to higher 
dimensions, as the famous Onsager solution of the 2D 
spin-1/2 Ising model shows [53|. Our approach has the 
advantage of offering a general formulation for any di- 
mension [53[. The problem within our approach is that 
the exact (in principle) solution is expressed in terms of 
a set of parameters not easily computable by means of 
the dynamics. The solution of this problem has been 
found for ID systems [53, Isl, [H, ]5^ and Bethe lattices 
[56|,|58|]; more complicated lattice structures are currently 
under investigation. Secondly, our approach provides us 
with an insight into the thermodynamic properties of the 
AL-EHM. Besides the exact determination of the zero- 
temperature phase diagram in the full range of the pa- 
rameters, one is also capable to analyze in detail the be- 
havior of several thermodynamic quantities such as the 
entropy, the specific heat and the charge and spin sus- 
ceptibilities. Furthermore, the study of these quantities 
allows one also to ascertain the different PTs occurring 
at zero temperature. Last but not least, the Green's 
function formalism and the equation of motion method 
allows one to exactly determine the energy spectra and, 
as a result, to identify the excitation energies contribut- 
ing to the specific heat. As a consequence, we are able 
to enlighten the issue of the nature of the peaks in the 
specific heat: we show that the double-peak structure is 
only due to charge excitations described by the Hubbard 
projection operators. 

The paper is organized as follows. In Sec. |TT1 upon 
introducing the Hubbard composite operators, we briefly 
review the analysis leading to the algebra closure and to 
analytical expressions of the retarded Green's functions 
(GF) and correlation functions (CF). Since the compos- 
ite operators do not satisfy a canonical algebra, the GF 
and CF depend on a set of internal parameters, leading 
only to exact relations among the CF. According to the 
scheme of the composite operator method [S^ , it is pos- 
sible to determine these parameters by means of algebra 
constraints fixing the representation of the GF. By fol- 
lowing this scheme, one can obtain extra equations clos- 
ing the set of relations and allowing for an exact and 
complete solution of the ID AL-EHM. In Sec. |TTT] we 
analyze the properties of the system at zero tempera- 
ture. We characterize the different phases emerging in 
the phase diagram drawn in the (U, n) plane, by study- 
ing the behavior, as functions of n and U, of the chemi- 
cal potential and of various local properties (double oc- 
cupancy, next-nearest correlation functions); the density 
of states and the charge and double occupancy correla- 
tion functions are also presented for various values of n 
and U. In the following, we shall denote with charge 
and double occupancy, the particle-particle and the dou- 



ble occupancy-double occupancy correlation functions, 
respectively. Section IIVI is devoted to the study of the 
finite temperature properties. Further to the study of all 
the quantities analyzed in Sec. IIIIl we also present results 
for the compressibility, the charge and spin susceptibili- 
ties, the specific heat, and the entropy. Finally, Sec. FVlis 
devoted to our conclusions and final remarks, while the 
appendices report some relevant but rather cumbersome 
computational details. 



II. EXACT SOLUTION OF THE ID ATOMIC 
EXTENDED HUBBARD MODEL 

The exact solution of the ID AL-EHM was presented 
in Ref . [54I , where the central issue was the case of van- 
ishing on-site potential {U — 0). In this section we shall 
briefly review the results previously obtained, providing 
a simpler framework of calculation and extending the so- 
lution to nonlocal correlation functions. This new frame- 
work allows one to study the properties of the system 
in the entire region of the external parameters n, T/V, 
U/V. 

By considering the atomic limit, the Hamiltonian of 
the extended Hubbard model in one dimension is given 
by 



H = Y^ [-fJ.n{i) + UD{i) + Vn{i)rf{i)\ 



(1) 



where n{i) — c^(i)c(i) is the charge density operator, 
c{i) [c^(i)] is the electron annihilation (creation) opera- 
tor - in the spinor notation - satisfying canonical anti- 
commutation relations. We use the Heisenberg picture: 
i = (i,t), where i stands for the lattice vector R^. The 
Bravais lattice is a linear chain of TV sites with lattice 
constant a. fi is the chemical potential; U and V are the 
strengths of the local and intersite interactions, respec- 
tively. D{i) — ri|(j)n|(i) = n{i)[n{i) — l]/2 is the double 
occupancy operator. Hereafter, for a generic operator 
$(«), we shall use the following notation 

<l>"(z,i) = ^ay<I'(j,i), 

3 

where aij is the projection operator over nearest neigh- 
boring sites 

"'J = ijJ2 e'^^^^'^'^aik) a{k) = cos(fca), 

k 

and the sum is over the first Brillouin zone. 



A. Eigenoperators and eigenvalues 

Upon introducing the Hubbard operators ^(i) — [n{i) — 
l]c(i) and rj{i) = n{i)c{i), one may define the composite 



field operator 



where 






(2) 



V'(«H^) 



V'^^HO^ 






\4'\^)j 



e(z)K(z)] 

/ r;(z) \ 
■q{i)[n°' [i)] 

Vr,WK(z)]4y 



(3) 



By exploiting the algebraic properties of the operators 
7i{i) and -D(j), it is easy to show that the fields 4''"^'^ (i) and 
ip^'^^i) are eigenoperators of the Hamiltonian ([T]) [52 |: 






(4) 



In Eq. (HI), e*^^) and e*^''^ are the 5x5 energy matrices 
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(5) 



£(") = 



\0 -3V fV -fV -fi+lOVj 

/U -fi 2V 

U-fi 2V 

U-n 2V 

U-fi 2V 

V -3V fy ~^v u-fi + iQV J 



whose eigenvalues, Em and Em\ are given by 



(6) 



E^^) = -^i + U+{m- 1)V 



{m = I,--- ,5}. 



(7) 



B. Retarded Green's functions and correlation 
functions 

The knowledge of a complete set of eigenoperators and 
eigenenergies of the Hamiltonian allows for an exact ex- 
pression of the retarded Green's function 

+ 00 



(2^) 



^^^-Mt-t')Q(s)^^^^ 



(8) 



and, consequently, of the correlation function 



+00 



(27r) 



dtje-'"^(*~*')c(^)(tj). 



(9) 



In the above equations, s = ^,77 and (• • ■ ) denotes the 
quantum-statistical average over the grand canonical en- 
semble. One finds 15211 



G^'^LO) 



5 



(s,m) 



„^1 Ul - E)n' + i5 



rn— 1 



(10) 



where T^n^ = 1 + tanh (^£;^ V2) , P = l/ksT, and the 
spectral density matrices cr^'*''") are given by 

^(.,m) ^5.(.)p(m)^ (11) 

where to = 1, • • • ,5; F^™^ are matrices of rank 5x5: 



.(1) 

Lk 



SLlSk. 



p(m) _ / TO - 1 
-■^ l,k — 



l+k-2 



.(s) 



where m — 2, ■ ■ ■ ,5 and the £,„ are given by 
e(^) = i (qI^^I - 25l['l + 35/(;j - 20l['l + 4li'^ 



5 ' 



e(^) = 



^ H ■ 



'1,3 ^ ^^1,4 ^-"i.s i ' 



S3'' = -e^i^s + 19/j;] ~ i6/j;i -t- 4/};^, 

■^4—0 \^^''l,2 ~ '-'1,3 + '^1,4 ~ ^-'1,5^ 

(12) 

-^a & ^'"^ ^^^ elements of the normalization matrix I^^' = 
({■0'-*''(i), '(/''''■' (i)}) which can be expressed as 

rW _ \(fe-l) 



A(^) = i(n(z)K(*)n. 



(13) 

Thus, the CFs depend on the external parameters n — 
{n{i)), T, U, V and on the internal parameters: /i, k^p\ 
and X'^pK It is easy to show that the CFs obey the fol- 
lowing self-consistent equations: 

C(«) + C^:^ - '^^'"'^ - At'^-i) (fc = l,-..5). (14) 

The number of these equations is not sufficient to de- 
termine all the internal parameters, and one needs other 
equations. This problem will be considered in the next 
subsection, where a self-consistent scheme, capable to 
compute the internal parameters, will be presented. 



C. Self-consistent equations 

The previous analysis shows that the complete solu- 
tion of the model requires the knowledge of the parame- 
ters /x, K^'^^ and A^'^'). These quantities may be computed 
by using algebra constraints and symmetry requirements 
|52 |. By recalling the projection nature of the Hubbard 
operators ^(i) and ry(i), it is easy to see that the following 
algebraic properties hold: 









(15) 



These are fundamental relations and constitute the basis 
to construct a self-consistent procedure to compute the 
various unknown parameters of the solution. One first 
fixes an arbitrary site, say i, of the chain, then splits the 
Hamiltonian JT]) in the sum of two terms: 



H,{i) = 2Vn{i)n°'{i), 



(16) 



and introduce the HQ{i) representation: the statistical 
average of any operator O can be expressed as 



(0) = 






(17) 



where (• • ■ )o,i stands for the trace with respect to 
the reduced Hamiltonian Ho{i)'. i.e., (■ ■ • )o,i — 



Tr{. 



.-l3Ho{'i.) 



}/Tr{6 



-iSHoii) 



}. In the following, we shall 



drop the index i since, by requiring translational invari- 
ance, all the sites are equivalent. As it is shown in 
Appendix |3 the parameters n^'^^ and A''"'-' can be writ- 
ten as a function of two parameters Xi = (n"(i))o and 
X2 — {D°'{i))o, in terms of which one may find a solution 
of the model. By exploiting the translational invariance 
along the chain, one can impose 



Hi)) 

{D{i)) 



(n"(z)), 
{D"{^)), 



(18) 



finding thus two equations allowing one to determine Xi 
and X2 as functions of fi. 

Xi ^ 2e''^(l -Xi- dX2){l + aXi + 0^X3) 

+ g/3(2M-c/) [2 + (d - l)Xi - 2dX2]{l + dXi + d^X2), 

X2 = p/(2a'-^) [1 + dXi - {2d + 1)^2] (1 + dXi + d^X2) 

- 2e^^{\ + d)X2{\ + aXi + a^Xa), 

(19) 



where a 



e-l3v _ 1 and d = e''^''^ 



1. The chemical 



potential /i can be determined by means of the equation 

I f X , - '/! X n m + a X 1 + «" X n ) 

(20) 



^^[{Xi-2X2){l + aXi 
+ 2X2(1 + ^X1+^2X2)], 



0^X2) 



where S = 1 - Xi + X2 + (Xi - 2X2) (1 + aXi + a^Xa) + 
X2{l + dXi + (PX2). Thus Eqs. ^ and ^ constitute 
a set of coupled equations allowing one to ascertain the 
three parameters /i, Xi, and X2 in terms of the external 
parameters of the model (n, t/, V , and T). Once these 
quantities are known, all the properties of the model can 
be computed. In particular, the local CFs can be ex- 
pressed in terms of the parameters D = {D{i)), k^p\ 
and \^P\ defined in Eq. ^, and tt^p) = {[D°'{i)]P}, 
J(p) = (n(i)[D"(i)]P)/2, e'-P^ = {D{i)[D"{i)]P)/2. In- 
deed, it is easy to derive the following relations for the 
two-point correlation functions 

{n{i)n{i + l)) = 2A(i) 

{n{i)n{i + 2)) = 2k(^) - n - 2D 

{n{i)D{i + l)) = 25(1) 

(n(^)i?(^ + 2)) = -«;(3) _ ^(2) + i^ ^^^^ 

o o 

{D{i)D{i + l)) = 2e''^^ 
{D{i)D{i + 2)) = 27r(2) -D, 
and for the three-point correlation functions 
{n{i)n{i + l)n{i + 2)) = 4A(2) - 2X^^'^ - 46'-^'^ 
{n{i)n{^ + l)D{r + 2)) = ^A^^) - 2A(2) + ^a'^) ^^^^ 
{D{i)n{i + l)D{i + 2)) ^ 4J(2) _ 2,5(1) 
{D{i)D{i + l)D{i + 2)) = 4e'(2) - 261(1^ 

The double occupancy, in terms of Xi and X2, is given 
by 



1 



D= :^[X2{l + dXi+d'^X2)]. 



(23) 



In Appendix [K\ we show how the parameters k^P', X^P', 
tt'^p\ S^p\ 9^p^ can be expressed in terms of the two basic 
parameters Xi and X2. At half filling, the system of 
equations PU)) and (P0)l can be solved analytically and 
one has 



D. Nonlocal correlation functions 

We shall now consider the problem of computing non- 
local correlation functions (Oi(«)02(j)), where Oi{i) and 
02{j) are generic operators and j is a site m steps 
{m > 1) away from the site i. To compute these quan- 
tities one needs only to know the local correlators ((2T|) 
thanks to some recurrence relations. We refer the in- 
terested reader to Appendix B for analytical details. In 
particular, for the charge and double occupancy correla- 
tion functions, one finds 



{n{i)n{j)) = n^ + Ap^' •''I 






(27) 



where the coefficients A, B, C, E, and the exponents p 
and q - all defined in Appendix B - are functions only of 
the short-range correlation functions X^^\ k,^^\ 9^^\ and 

^(2). 



E. Density of states 

By noting that the cross GFs {■ip^^\i,t)->p^'^^\i,t)) van- 
ish, the electronic density of states (DOS) is given by 



N{u;)=(--)lm[G[^^{co) + G[\\u 



J2 [^li^sico - £;(«)) + E(„^)5(^ - £;i^))]. 



(28) 



As one can immediately see, the DOS is a sum of delta's 

(s) 

functions weighted by the corresponding Sm s. It is 

(s) 

worthwhile to observe that the expressions of the Sm s 
(fT2)) lead to the following sum rule: 



EPi'^ + 4'"] 



dLuN{uj) = 1. (29) 



^l = -u + 2V, 



- v/(l + 2Gi^ + if2)2 _ 8GK{K - 1)2]}, 
X2 = 7^ rrrjl + 2GK + K^ 

- \/(l + 2Gif + K^f - 8GK{K - 1)2}, 



(24) 



(25) 



where K ~ e ^^ and G ~ e^^'"^. As a consequence, one 
finds 



D = 



1 



2 + 2G(1 - a2A:2 



2 V„^2 



(1) _ 1 - 2a^X2 + {l-AK + 2^2 ^ j^4)^2 (26) 

" 2[1- 202^2 +a2(l + if2)^2] ■ 



In the following, we shall use a Lorentzian broadening for 
the delta functions S{uj — ujo) = e/[{uj — cjo)^ + e^l^r with 

£ = 0.1. 



F. Summary 

Summarizing, in this section we have shown that the 
model can be exactly solved by using the Green's function 
formalism and the equation of motion method. The cen- 
tral point in this approach lies in the fact that there exists 
a closed set of eigenopcrators of the Hamiltonian ^ [see 
Eqs. ([3])] allowing for an exact determination of the cor- 
relation functions, as it is shown in Sees. IIIBl IIIDI and 
in the Appendices. Exact expressions for the GFs and 
CFs can be written in terms of few local correlators. By 
using algebraic properties of the relevant operators and 
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FIG. 2: The charge (a) and the double occupancy (b) corre- 
lation functions as functions of the distance m = \i ~ j\ for 
t/ < at T = 0, V = 1 and n=0.25, 0.5, 0.75, 1. 



FIG. 1: The phase diagram in the plane {U,n) at T = and 
V = 1. One can identify four phases. 



symmetry properties, these correlators can be expressed 
in terms of the chemical potential /i and of the parame- 
ters Xi and X2. These three quantities are determined 
as functions of the external parameters n, [/, V^ and T 
by solving the system of equations ((T9)) and (|20|) . Once 
these quantities are known, all properties of the model 
can be computed. In the next section we shall present a 
comprehensive study of the properties of the system at 
zero-temperature . 



III. ZERO TEMPERATURE RESULTS 

In this section we derive the phase diagram of the AL- 
EHM in the (J7, n) plane. More specifically, we shall nu- 
merically solve the set of equations ([T51) and (PO)) at T = 
and study the behavior of relevant correlation functions 
and of the density of states in order to envisage the distri- 
bution of the particles for different densities. The results 
obtained are displayed in Fig. [TJ one may distinguish 
four different phases, according to the values of the in- 
ternal parameters. 

Phase (a). The phase (a) is observed in the region 
(0 < n < 1, f7 < 0) and is characterized by the following 
values of the parameters: 



Xi = 


2n 


2-n 


x^ 


n 


-^2 


2-n 


D^ 


--n/2, 


J2) _ 


2n 



_(3) _ 2n(l -t- n) A = {2-n)n, 

2-n ' B = E = 0, 
^(,)^Ml + 3n), ^^{2-n)n 



r(2) 



2-n 
n 

2(2 -n)' 
2-n A^*^) = e'^i' = 0, 



4 
n 



P = 



2-n 



q = 0. 



The chemical potential takes the value U/2 for n < 1, 
whereas at n = 1 one naturally finds fi = 2V + U/2. The 
attractive on-site potential favors the formation of doubly 
occupied sites. At the same time, the repulsive intersite 
potential disfavors the occupation of neighboring sites. 
The expectation values D = n/2 and (n(i)[n"(i)]'^) = 



suggest that, in the thermodynamic limit, there are nei- 
ther singly occupied sites nor neighboring sites occupied. 
For n < 1 one does not find any charge ordered state and 
the energy necessary to add two electrons is 2^ = U. At 
half filling (n = 1), the chemical potential jumps from 
the value fi — U/2 to the value fj, = 2V + U/2, as re- 
quired by particle-hole symmetry. Indeed, one needs to 
provide an energy amount of 21^ for an extra electron to 
occupy a nearest neighbor site. At n = 1 half of the sites 
are doubly occupied: by varying the particle density, the 
system undergoes a FT to an ordered state characterized 
by a checkerboard distribution of doubly occupied sites. 
This can be clearly seen in Figs. [2^-b, where we plot the 
correlation functions 

{n{i)n{j)) =n^ + n{2 - n)(-l)™ e~™/«", 

{D{^)D{J)) =D' + !li^(-l)"e-'"/«. ^^^^ 

as a function of the distance m — \i — j\. At half filling, 
both CFs present a two-site periodicity. The correlation 
length ^„ in Eqs. (|5D|) is defined as 



^n 



In 



(31) 



When n < I the parameter \p\ is less than 1. Thus, the 
CFs (pO)) decrease exponentially, oscillating around the 
ergodic values n^ and D^, respectively. On the other 
hand, at half filling \p\ = 1: the correlation length di- 
verges, a long-range order is established, and both CFs 
present a two-site periodicity. 

The distribution of the particles in the phase (a) is 
shown in Fig. [3l where we report just one possible con- 
figuration. The configurations in the ground state will 
have no doubly occupied sites next to each other. When 



FIG. 3: Distribution of the particles along the chain by vary- 
ing the particle density at T = and U < 0. White and black 
squares denote empty and doubly occupied sites, respectively. 
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FIG. 4: (a) The density of states at zero temperature, V — 1 
and U = —1.5 is reported as a function of the frequency for 
n = 0.25, 0.5 and 0.75. (b) The density of states at zero 
temperature, V = 1 and n = 1 is reported as a function of 
the frequency for U = —1,-2 and —3. 
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FIG. 5: The charge correlation function {n{i)n{j)) for [/ > 
is reported as a function of the distance m = \i — j\ at T = 0, 
V = 1 and n = 0.3, 0.4 and 0.5. 



n < 1, there is no ordered pattern in the distribution of 
the particles, whereas, for n = 1, one observes the well- 
known checkerboard distribution of doubly occupied sites 
[la |. In Figs. [4^-b, we plot the density of states at zero 
temperature and attractive U, for n < 1 and n = 1, re- 
spectively. For n < 1, the DOS exhibits four peaks at 



pM 



AO 



(«) 



energies £['" = U/2, E[^' = -U/2, E^^^' = -U/2 + 2V, 

E^^' = — C//2 + AV: the electrons can occupy only states 
with these associated energies. For n — 1, the excita- 
tion spectrum contains only two allowed values, leading 
to a DOS with two peaks at e['^^ = U/2 - 2V, and 

-Eg — -U/2 + 2V. By varying the particle density, 
the weights of the peaks change as 



1 2 



Ei' 



(j)_2(l-n) 



.(C) 



2n(l - n) 



2~n 

2 



5 2(2- 



The weights relative to the other energies are zero: E2 

0. 



(«) 



v(C) — v^**) — v('') — v^'') — v^")) 

Z-IA — Zjo — Zjq — Zj^ — Zj 



Phase (h). The phase (6) is observed in the region (0 < 
n < 0.5, [/ > 0) and is characterized by the following 
values of the parameters 



Xi = - , 

1 — n 

^(2) = 0(1) = 
A« = 0, 

and 



c(2) 



^(3) = 



0, 



c(4) 



2(1 -n)' 
n(l + 2n) 


4(1 - n) 



8(1 



A 
B 

P 

q 



n(l — n), 

E = C^O, 
n 



1-n 



0. 



{n< 0.5) 

U {n = 0.5 and < C/ < 2V) 

2V (n = 0.5andC/> 2V). 



The repulsion between electrons on the same site and 
on neighboring sites, leads to a scenario where the elec- 
trons will tend to singly occupy non- neighbor sites, 
as confirmed by the expectation values D = and 



{n{i)[n°' {{)]'') = 0. For particle densities less than quar- 
ter filling, i.e., n < 0.5, there is no cost in energy to add 
one electron, thus 11 — 0. At quarter filling, one has to 
distinguish two cases: (i) when < U < 2V, the chem- 
ical potential jumps from the value fi ~ to the value 
ji = U\ (ii) for U > 2V, the chemical potential presents 
a jump from /x = to /i = 2V. As a result, for n > 0.5, 
as illustrated in Fig. [TJ there are two distinct phases, (c) 
and (d). When < U < 2V, by adding one electron the 
system enters phase (c) , characterized by the presence of 
both singly and doubly occupied sites. The minimization 
of the energy leads to a cost in energy oi f.i = U to ac- 
commodate an extra electron. On the other hand, when 
U > 2V, by adding one electron, the system is driven 
into phase (d), characterized by the presence of neigh- 
bor arbitrary-spin singly occupied sites, requiring thus 
an energy cost of ^ — 2V. 

Right at quarter filling, half of the sites are singly 
occupied: by varying the particle density, the system 
undergoes a FT to an ordered state characterized by 
a checkerboard distribution of arbitrary-spin singly oc- 
cupied sites. This is clearly shown in Fig. [51 where 
we plot the charge correlation function {n{i)n{j)) — 
m? + n(l — n)(— 1)™ e"™/'^" as a function of the rela- 
tive distance: when n < 0.5 one finds that the parame- 
ter \p\ appearing in Eq. (P7)) is less than 1. As a result, 
{n{i)n{j)) decreases exponentially, oscillating around the 
ergodic value n^. At quarter filling one finds |p| = 1: the 
correlation length pip diverges and a long-range order 
is established, as evidenced by the two-site periodicity of 
{n{i)n{j)). 

The DOS at zero temperature and repulsive U is shown 
in Fig. Hi for n < 0.5. In Fig. [Sb we plot the DOS 
at n = 0.5 for U < 2V (continuous line) and U > 2V 
(dotted line). For n < 0.5, the DOS exhibits four peaks 
at the energies e['^^ = U, e[^^ 



= 0,E. 



(«) 



V, eP = 2V. 



For n = 0.5 and < U < 2V, the electrons can occupy 



only three states at the energies El' — 0, E^ 

E^ = -U + 2V. For n = 0.5 and U > 2V, the DOS 

still exhibits three peaks but now at the energies El'' — 

U-2V,e[ 



(€) 



-u, 



2V, E^ = 0. By varying n the weights 
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FIG. 6: The density of states N{ll}) as a function of the fre- 
quency ll! at zero temperature, V = 1 and (a) U — 1/2 for 
n=0.2, 0.3 and 0.4; (b) n = 0.5, for U = 3/2 and 3. 



FIG. 7: The charge (a) and the double occupancy (b) cor- 
relation functions as functions of the distance m = \i — j\ at 
T = 0, F = 1 and n =0.5, 0.75 and 1. 



of the peaks change as 



Av) _ 
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2n(2n- 1) 
1-n 



(5) _ 2-7n(l-w) 
1 " 2(1 -n) 
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1 — n 



The weights relative to the other energies are zero. 

Phase (c) The phase (c) is observed in the region (0.5 < 
n < 1, < U < 2V) and is characterized by the following 
values of the parameters: 



Xi = 2n, 
X2^2n- 1, 
D^n-1/2, 

7r(2) ^n{n- 1/2), 
A(fe) = 0(1) = 0, 



.(2) - 



1(271^ + 3n-l), 



t(3) = (18^2 + n- 3), 



c(4) = i 



(110n2-45n-l) 



and 



B = E = -2n'^ + 3n-l, ^ 



C 



\I2)\ 



The chemical potential takes the value U for 0.5 < n < 1 
and 2V -I- {//2 at n = 1. The repulsive on-site and in- 
tersite potentials force the electrons to prefer seating at 
different sites and disfavor the occupation of neighboring 
sites. However, when the intersite interaction dominates 
the on-site one, the minimization of the energy requires 
the extra electrons above n = 0.5 to occupy nonempty 
sites. This situation is well described by the expectation 
values D = n - 1/2 and {n{i)\n°' {i)f) = 0. Thus, in 
the thermodynamic limit, there can be singly and dou- 
bly occupied sites and no neighboring sites are filled. By 
increasing n, the number of doubly occupied sites in- 
creases linearly. In this region, the correlation length 
(PIj) diverges: one finds a charge ordered state character- 
ized by a checkerboard distribution of alternating empty 
and occupied (either singly or doubly) sites. For n < 1 
the energy necessary to add one electron is fj, = U. Upon 
increasing the particle density, at half filling the checker- 
board distribution is realized with only doubly occupied 



sites alternating with empty ones and, correspondingly, 
the chemical potential jumps from the value fj, = U to 
the value fi — 2V + U/2. Since all the nonempty sites 
are already doubly occupied, an extra electron would go 
to occupy an empty site with a cost in energy of 2V. 
The existence of the ordered states is endorsed by the 
behavior of the charge and double occupancy CFs. 



{n{i)n{i + m)) =n2[l + (-l)'^ 
{D{i)D{i + m)) ^ D^[l + i-iy 



TO > 1 



and reported in Figs. [7^ and [7)3. One clearly observes a 
two-site periodicity of (n(j)n(i-l-m)) and {D{i)D{i + m)), 
leading to an ordered charge distribution: one every two 
sites is at least singly occupied. The double periodicity 
of the double occupancy CF is, of course, to be under- 
stood as an average over all the possible distributions of 
doubly occupied sites that are, of course, energetically 
equivalent. 

The distribution of the particles in the region (0 < 
n < 1, < U < 2V) is drawn in Fig. [HJ where we 
represent just one possible configuration. At quarter fill- 
ing the ground states of the Hamiltonian are checker- 
board configurations where empty sites alternate with 
sites occupied by one particle of arbitrary spin. By in- 
creasing n, the charge order persists and one observes 
an increase of the number of doubly occupied sites, as 
illustrated in Fig. [51 Finally, at half filling, one ob- 
serves a checkerboard distribution of doubly occupied 
sites [2J]. The density of states at zero temperature is 
plotted in Fig. [9^. for 0.5 < n < 1 and U = V/2 and 
in Fig. [Hb for n = 1 and different values of U. One 
observes that, when n <1, the DOS exhibits five peaks 
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at the energies El ' — 0, E^ 



Ai) 



-u, eP 



-U + 2V, 



FIG. 8: Distribution of the particles along the chain by vary- 
ing the particle density at T = and < U < 2V. White 
and black squares denote empty and doubly occupied sites, 
respectively. Grey squares denote arbitrary-spin singly occu- 
pied sites. 
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FIG. 9: The density of states N{ll}) as a function of the fre- 
quency u at zero temperature, V = 1 and (a) [/ — 1/2 for 
n = 0.6, 0.7, 0.8, 0.9; (b) n = 1, for [/ = 1/2, 1 and 3/2. 



FIG. 11: The density of states N{lu) as a function of the 
frequency to at zero temperature, V = 1 and (a) U — 2.1 
for n = 0.6, 0.7, 0.8, 0.9; (b) n = 1, for [/ = 5/2, 7/2 and 
t/ = 9/2. 
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The weights relative to the other energies are zero. 

Phase (d) The phase (d) is observed in the region 
(0.5 < n < 1, [/ > 2V) and is characterized by the 
following values of the parameters, 
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In this region the on-site interaction dominates the in- 
tersite one. The minimization of the energy requires the 
electrons not to be paired and allows for the occupation 
of neighboring sites: D = and (n(z)[n"(i)]'') ^ 0. For 
n < 1 the energy necessary to add one electron is /i = 2V. 
At half filling the chemical potential jumps from the value 
fi — 2V to fj, — 2V + U/2. This is because an extra elec- 
tron would go to occupy a singly occupied site, requiring 
thus an energy equal to U. The behavior of the charge 
correlation function is reported in Fig. (TU] as a function 
of the relative distance. One observes that increasing 
the particle density, {n(i)n{j)) becomes more and more 
uniform losing the two-site periodicity peculiar of the re- 
gion < [/ < 2V . At half filling the system presents an 
homogeneous distribution, with all singly occupied sites 

The density of states at zero temperature and U > 2V 
is plotted in Fig. [TTk for 0.5 < n < 1 and in Fig. [Tib 
for n = I and different values oi U {U > 2V). For n <1, 

the DOS exhibits five peaks at the energies E^ ' — —2V, 

For 



e: 



(?) 



-V, E. 



(?) 



0, E. 



iv) 



u -V, e: 



iv) 



u. 



1, the electrons can occupy two states at the en- 



ergies E. 



(?) _ 



;.('?) _ 



= -U/2, E^'" = u/2. 



The distribution of 



the particles in the region {0 < n < I, U > 2V) is pre- 
sented in Fig. [I3 where just one possible configuration 



n=0.2S 
tt=0.5 

}i=0.75 



FIG. 10: The charge correlation function {n{i)n{j)) for U > 
2V as a function of the distance m = |j — j| at T = 0, V = 1 
and n = 0.6, 0.8 and 1. 



FIG. 12: Distribution of the particles along the chain by vary- 
ing the particle density at T = and U > 2V. White and 
grey squares denote empty and arbitrary spin singly occupied 
sites, respectively. 
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(a) (b) 

FIG. 13: (a) The chemical potential fi as a, function of the 
particle density for V = 1, T = and U — —1.5, 1.5 and 
2.5. (b) The chemical potential as a function of the on-site 
potential U for V — 1, T = in the region < n < 1. 







(a) (b) 

FIG. 14; (a) The double occupancy as a function of the par- 
ticle density for V = 1, T = 0, (/ < (dot-dashed line), 
< U <2V (dotted line) and U > 2V (continuous Une). (b) 
The double occupancy as a function of the on-site potential 
U for V — 1, T — and various values of n. 



is shown. As one can clearly see, in region (d) there are 
no charge ordered states but at ti = 0.5. 

One can now investigate some thermodynamic quan- 
tities, whose behavior can also be relevant to ascertain 
some of the results previously obtained. In particular, we 
shall investigate the behavior of the chemical potential 
and of the double occupancy as a function of the particle 
density and of the on-site potential. In the former case, 
we shall consider the full range < n < 2 to better rec- 
ognize the peculiar behaviors of /i and D. In Fig. [T5b we 
plot the chemical potential as a function of n for differ- 
ent values of U: fi exhibits plateaus and discontinuities. 
In particular, (i) for [/ < 0, /i takes the constant value 
/i = C//2 in the region < n < 1, corresponding to phase 
(a). At n = 1 there is a discontinuity corresponding 
to the PT to a checkerboard distribution of doubly occu- 
pied sites; (ii) for < U < 2V, fi takes the constant value 
/^ = in the range < n < 0.5, corresponding to phase 
(6). At n = 0.5 there is a discontinuity corresponding 
to the phase transition to the checkerboard distribution 
of singly occupied sites. In the range 0.5 < n < 1, /i 
takes the constant value jjl = U , corresponding to phase 
(c). At n = 1 there is a discontinuity corresponding to 
the PT to the checkerboard distribution of doubly occu- 
pied sites; (iii) for U > 2V, jjl takes the constant value 
jj, = 2V in the region 0.5 < n < 1, corresponding to 
phase (d). At n = 1 there is a discontinuity correspond- 
ing to the PT to the homogeneous distribution of singly 
occupied sites. When plotted as a function of U [see Fig. 
[T3b]. the chemical potential shows a linear behavior for 
n = 1, as required by the particle- hole symmetry. For 
all other values of the filling, ^ presents a linear behavior 
for attractive on-site interaction. Then, for positive U, 
fi presents a plateau (^ = 0) when n < 0.5, whereas, for 
0.5 < n < 1, it increases linearly with U up to U — 2V, 
where it takes the constant value fi = 2V. As it is evident 
in Fig. [Hb . the double occupancy shows three different 
behaviors according to the value of U. For < n < 1, the 
double occupancy is always zero when U > 2V, it shows 
a linear behavior for negative values of U and presents 
a discontinuity at quarter filling, when < U < 2V. In 
Fig. [Mb we plot the double occupancy as a function of 
the on-site potential U for several values of n: one ob- 
serves a one-step or a two-step behavior depending on 
the particle density. 



In closing this section, it is worthwhile to mention that 
at zero temperature all the phases (with the exception 
of the case n = 1 and U < 2V) exhibit a macroscopic 
degeneracy growing exponentially with the volume of the 
lattice. 



IV. FINITE TEMPERATURE RESULTS 

In this section we shall investigate the finite temper- 
ature properties of the AL-EHM. We study the behav- 
ior of the system as a function of the parameters n, T, 
and U and again we take 1^ = 1 as the unit of energy. 
We shall focus our attention on several thermodynamic 
quantities: namely, the chemical potential /x, the dou- 
ble occupancy D, the nearest-neighbor charge correlation 
X°' = 2A(^) = {n{i)n"{i)}, the thermal compressibility 
kC^) — (1/n.^) • dn/dii, the density of states, the charge 
{n{i)n{j)) and double occupancy {D{i)D{j)) correlation 
functions, the charge Xc and the spin x^ susceptibilities, 
the internal energy E = UD + Vx"', the specific heat 
C = dE/dT, and the entropy S. Under the transforma- 
tion n — > 2 — ri, these quantities transform as: 



^i{2- 


- n) 


= W + U -^l{n) 


D{2- 


- n) 


= l-n + D{n) 


X"(2- 


- n) 


^A{l-n) + x"{n) 


Xc(2- 


- n) 


= Xc{n) 


Xs{2- 


- n) 


= Xs{n) 


E{2- 


- n) 


= (U + AV){1 - n) + E{n) 


C(2- 


- n) 


= C{n). 



(32) 



Therefore, unless otherwise stated, we shall limit the 
analysis to the range < n < 1. In order to charac- 
terize the various features of the thermodynamic quanti- 
ties, in the following we shall distinguish three intervals 
along the U line (at the borders of which zero temper- 
ature transitions occur): namely, U < 0, < U < 2V, 
and U > 2V. 
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(c) 



FIG. 15: The chemical potential /i as a function of the particle 
density n for various temperatures and (a) U < 0, (b) < 
U < 2V, (c) U > 2V. 
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FIG. 16: The chemical potential /i as a function of the on-site 
interaction f/ for n = 0.25, 0.5, 0.75 and for T = 0.2 (solid 
line) and T = 0.5 (dotted line). 
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The behavior of the chemical potential as a function 
of the particle density n is reported in Figs. [TSk-c as the 
temperature varies, for the three regions C/ < 0, < C/ < 
2V and U > 2V, respectively. One can immediately see 
that /x is always an increasing function of n, hinting at 
a thermodynamically stable system for all ranges of the 
parameters n, T, and U. We consider the full range < 
n < 2 to better appreciate the jumps in fi as the particle 
density is varied. An interesting feature is the presence 
of crossing points in the chemical potential curves, when 
plotted for different temperatures and U < 2V. More 
precisely, one observes crossing points at the values n = 
0.5, n = 1, and n = 1.5, where PTs are observed at 
T — 0. In between these values, for temperatures close 
to zero, the chemical potential presents plateaus. In Sec. 
II V HI we shall see that nearly universal crossing points are 
also exhibited by the specific heat. The value n — 0.5 is 
a turning point for the derivative dji/dT: for n < 0.5 
(n > 0.5), /i is a decreasing (increasing) function of T. 
In Fig. \W\ we show the chemical potential as a function 
of the on-site interaction U ior n — 0.25, 0.5, 0.75 and 
for two temperatures. For negative U, fi increases by 
increasing U following approximately a linear law (/i = 
[7/2 for T — 0). When U > 0, ^ deviates from the hnear 
behavior and, for U « 2V, tends to a constant value 
depending both on n and T. The latter behavior has 
been observed also for the 2D Hubbard model (see Ref. 
[60| and references therein). 



(c) 



FIG. 17: The double occupancy as a function of the tem- 
perature for different values of the particle density and (a) 
U = -1/2, (b) U^l,{c)U = 5/2. 



B. The double occupancy 

In Figs. [T7k-c the temperature dependence of the dou- 
ble occupancy D is shown for the three regions U < 0, 
< U < 2V and U > 2V, respectively. In the first 
region (attractive U), D is always a decreasing function 
of the temperature since the thermal excitations tend to 
break the doublons, i.e., the charge carriers of doubly oc- 
cupied sites. In the second region {0 < U < 2V), the 
double occupancy is an increasing (decreasing) function 
for n < 0.5 (n > 0.5). This is easily understandable if 
one thinks of the zero-temperature scenario: for filling 
less than one-quarter, the thermal excitations may favor 
the formation of doublons, against the repulsive interac- 
tion U. For filling greater than one-quarter, the thermal 
excitations may favor the formation of neighboring singly 
occupied sites, against the repulsive interaction V. In the 
third region {U > 2V), the double occupancy is always 
an increasing function of T, since the thermal excitations 
favor the formation of doublons, against the repulsive in- 
teraction U. A similar behavior of the double occupancy 
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FIG. 18: The double occupancy as a function of the on-site 
interaction U for n — 0.25, 0.5, 0.75, 1 and for two different 
temperatures. 



by varying the temperature has been found in a strong 
couphng analysis of the EHM at half fiUing [43] • 

It is interesting to plot the double occupancy as a func- 
tion of the on-site potential U. As expected, at all tem- 
peratures, D is always a decreasing function of U. As 
one can see from Fig. [THl the discontinuities observed 
at zero temperature for U — and U = 2V are more 
and more smoothed as the temperature is increased. In 
the limit C/ — > — oo the double occupancy tends to the 
constant value D = n/2, independently of the tempera- 
ture: all the electrons are paired. In the opposite limit, 
namely [/ — > oo, the on-site repulsion prevents the double 
occupations of the sites. 



The nearest-neighbor charge correlation 
function 



As discussed in Sec. IIIIl there are two interesting quan- 
tities which mainly characterizes the state of the system: 
one is the double occupancy, discussed in the previous 
subsection, the other is the nearest-neighbor charge cor- 
relation function x" = {n{i)n°'{i)) = 2A'^'^^. The behav- 
ior of x" as a function of the temperature is shown in 
Figs. [Tgk-c. for the three regions U < 0, < U < 2V 
and U > 2V, respectively. As one would have expected, 
the behavior of x" is opposite to the one of D in some 
regions of the parameters: for instance, for attractive 
on-site interactions, thermal excitations tend to break 
the doublons and to consequently favor the occupation 
of neighboring sites. For U < 2V, x" is always an in- 
creasing function of T and vanishes in the limit T — > 0, 
where occupation of neighbor sites is clearly disfavored. 
As one can see in Figs. [T5k-c. the thermal excitations 
may favor the occupation of neighboring sites and, in the 
limit T — > oo, x" tends to the value n^ for all values 
of U and n. When the on-site potential is greater than 
2V, the behavior of x" depends on the filling: it is an 
increasing (decreasing) function for n < 0.5 {n > 0.5) at 
low temperatures. This different behavior can be under- 
stood if one looks at Fig. [T^ it is clear that for lower 
fillings thermal excitations are more effective in augment- 
ing the nearest-neighbor charge correlations with respect 
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(c) 



FIG. 19; The nearest-neighbor charge correlation function 
X™ ~ {n{i)n°' (i)) as a function of the temperature for different 
values of the particle density and (a) U = —1/2, (b) U — 1, 
(c) U = 5/2. 



to higher fillings, where electrons are already clustered 
together. 

In Fig. [50] we plot x" ^s a function of the on-site 
potential: one observes that x" is always an increasing 
function of U , and tends to a constant value (which de- 
pends on n), independently of the temperature, in the 
limit of C/ — > -|-oo. 
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FIG. 20: The nearest-neighbor charge correlation function 
X" ~ {n{i)n°' (i)) as a function of the on-site interaction U 
for n = 0.25, 0.5, 0.75, 1 and for T = 0.2 (solid lines) and 
r = 0.5 (dotted lines). 



D. The compressibility 

In this subsection we study the thermodynamic com- 
pressibility, defined as k = {dn / d fi) / -n? . In Figs. [5Tb-c 
we plot K as a function of the filling n, for three repre- 
sentative values oiU {U ^ -V, U ^V imdU ^ 2.W, 
respectively) and for different temperatures. To avoid di- 
vergences in the limit n — > 0, we plot only the derivative 
(dn/dfj.). In the first region {U < 0, phase (a)), the com- 
pressibility increases by increasing n up to a maximum 
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FIG. 21: dn/dn as a function of the particle density n for 
various temperatures and (a) f/ < 0, (b) < U < 2V, (c) 
U>2V. 



in correspondence of quarter-filling; then it decreases ex- 
hibiting a minimum at half-filling, where a charge ordered 
state is established at T = (see Fig. [3]). In the second 
region {0 < U < 2V, phases (b) and (c)), k has a dou- 
ble peak structure with two maxima at n = 0.25 and 
n = 0.75 and two minima at n = 0.5 and n = 1; the 
latter points correspond to the situation where a charge 
ordered state (checkerboard distribution for singly and 
doubly occupied sites, respectively) is observed at zero 
temperature (see Fig. [5|). In the third region {U > 2V, 
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FIG. 23: The charge (a), (b), (c) and the double occupancy 
(d), (e), (f) correlation functions as a function of the distance 
m — \i ^ j\ for pertinent values of the temperature, V = 1 
and U = 1/2. The value of the filling grows from left to right. 



the filling, there are only singly occupied sites (see Fig. 
[T^ . In Figs. [22^-c we plot {dn/dfj,) as a function of the 
temperature, again in the three regions oi U {U = —V, 
U = V and U — 2.1V, respectively) and for different val- 
ues of the filling {n = 0.25, 0.5, 0.75, 1). For all values of 
U and n, in the limit of large temperatures, the compress- 
ibility vanishes following the law n ex T^^, but becomes 
strongly filling dependent as T is reduced. In particular, 
in the limit T — > and in the regions where no CO state 
is observed, n increases by lowering T diverging at T = 0. 
On the other hand, in the regions where a charge ordered 
state is present, the compressibility exhibits a maximum 
at a certain temperature, then, further decreasing T, k 
vanishes. 



(c) 



FIG. 22: dn/dn as a function of the temperature for various 
fillings and (a) C/ < 0, (b) < C/ < 2V , (c) U > 2V. 

phases (b) and (d)), n has a behavior similar to the one 
exhibited in the second region, with the difference that 
the two peaks have now the same height. This difference 
is due to the fact that in the phase (d), for all values of 



E. Correlation functions at finite temperature 

At low but finite temperatures, the charge and 
double occupancy correlation functions {n{i)n{j)) and 
{D{i)D{j)) still show a finite-range order evidenced by 
a two-site periodicity. At quarter filling, one clearly ob- 
serves in Figs. [23k and [23H that one every two sites is 
singly occupied: the two-site periodicity of the charge CF 
decreases with increasing T and distance m. The double 
occupancy CF is still very close to zero hinting at the ab- 
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sence of doubly occupied sites also at low temperatures. 
Upon increasing the filling, for instance at n — 0.75 (see 
Figs. [23b andESfe), both {n{i)n{j)) and {D{i)D{j)) ba- 
sically retain the same zcro-tcmperature behavior also at 
low temperatures. At half filling the double periodicity 
of the CFs persists even at higher temperatures, as evi- 
denced in Figs. [^5b and[ 



F. Charge and spin susceptibilities 

Another important quantity useful for studying the 
critical behavior of the system is the susceptibility. In 
this subsection we shall compute the charge and spin sus- 
ceptibilities. 



1. Charge susceptibility 

Let A(/c) be the Fourier transform of the charge corre- 
lation function 



(n(j>(0) = 



27r 



r/a 



— ■K I a 



dke'''^'-^^A{k). 



From Eq. ([27]) it is easy to derive 



2^ 2e„„^ , A{l-p^) 



A(k) = —n'5(k) + 

^ ' a 1+p^ 



2pcos(fca) 



5(1 -g2) 



1 + g^ — 2qcos{ka) 



Then, it immediately follows that the charge static sus- 
ceptibility Xc is given by 



Xc = A(0) = Nn^ 



A{l+p) Bjl + q) 
(l-p) ^ (1-g) ' 



(33) 



where N is the number of sites. It is worthwhile to recall 
that Xc can be also computed by means of the thermo- 
dynamics through the formula 



Xc 



Nn^ 



]_dn 



(34) 



The difference between the above two expressions for Xc 
lies in the fact that ([55)1 requires the knowledge of a two- 
point electronic Green's function, while (j34[) requires only 
the knowledge of a one-point electronic Green's func- 
tion. In condensed matter theory there are many ex- 
amples where a physical quantity can be computed in 
different ways by using GFs of different hierarchy (see, 
for example, Sec. 1.5 in Ref. [5^). Of course the re- 
sult should not depend on which path of calculation one 
has followed. This is true if one is capable to exactly 
solve the model. However, in many practical situations 
this is not possible and one needs to use some approx- 
imations. As a consequence, one may obtain different 
results, because approximation methods usually perform 
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FIG. 24: The charge susceptibility as a function of the particle 
density for V = 1 and different values of the temperature, and 
for (a) [/ = -1 ; (b) t/ = 1; (c) U = 2.1. 



in different ways according to the rank of the GFs under 
analysis. Actually, by comparing different expressions of 
the same quantity one may envisage the effectiveness of 
an approximation. For the present model, it is not diffi- 
cult to show that Eqs. ([55]) and (|34p give the same result 
since the model has been exactly solved. In the following, 
we will use one or the other according the convenience of 
numerical calculations. At T = 0, the charge susceptibil- 
ity can be analytically computed since one knows all the 
coefficients A, B, p, and q, as they are given in Appendix 
mil In particular, one can compute Xc in each region of 
the phase diagram: 



Xi'^'=n(l-n)(2-n), 
X?)=n(l-n)(l-2n), 



.(<=) 



(1 -n)(2n- 1), 



Xi''^=n{l-n)i2n~l). 

The superscript in the above expressions labels the 
different phases at T = 0. In Figs. [24k- c we plot Xc 
as a function of the filling n for three fixed values of 
U {U ^ -V, U = V andU = 2AV, respectively) and 
for different temperatures. One immediately sees that, at 
T = 0, the charge susceptibility is finite except at n = 0.5 
and n = 1, where, at fixed on-site potential, PTs occur. 
In the first region [U < 0, phase (a)], the charge suscep- 
tibility increases by increasing n and has a maximum in 
correspondence of quarter filling. Further increasing n, 
Xc decreases and vanishes at half filling, where a charge 
ordered state is estabfished at T = (see Fig. [3]). In 
the second region [0 < U < 2V, phases (6) and (c)], for 
low temperatures, Xc has a double peak structure with 
two maxima around n ~ 0.25 and n = 0.75 and two 
minima around n = 0.5 and n = 1; the latter points 
correspond to the situation where a charge ordered state 
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FIG. 25: The charge susceptibihty as a function of the tem- 
perature for 1/ = 1, n = 0.25, 0.5, 0.75, 1 and for (a) U = -\ 
; (b) [/ = 1; (c) V = 2.1. 



FIG. 26: The spin susceptibihty as a function of the temper- 
ature for 1/ = 1, n = 0.25, 0.5, 0.75, 1 and for (a) [/ = -1 ; 
(b) (7 = 1; (c) U = 1.1. 



(checkerboard distribution of singly and doubly occupied 
sites, respectively) is observed at zero temperature (see 
Fig. [5]). In the third region \1V < U, phases (&) and 
(d)], Xc shows a behavior similar to the one exhibited 
in the second region, the difference being that the two 
peaks have the same height. This difference is due to the 
fact that in the phase [d) one has only singly occupied 
sites for all values of the filling (see Fig. [T^ . For high 
temperatures, Xc is always an increasing function of n. 
In Figs. [^5b-c we plot Xc as a function of the tempera- 
ture, again in the three regions oiU {U = —V, U = V 
and U = 2.iy, respectively) and for different values of 
the filling {n = 0.25, 0.5, 0.75, 1). For aU values of U 
and for large temperatures, the charge susceptibility in- 
creases and tends to a constant value in the limit of T 
infinite. At low temperatures, one observes a different 
behavior: at T = 0, Xc is finite for values of the particle 
density corresponding to nonordered states. Further in- 
creasing T, Xc decreases, it reaches a minimum and then 
increases again. On the other hand, Xc decreases with 
T and vanishes at T == for fillings corresponding to 
charge ordered states. In the linfit of high temperatures, 
the charge susceptibility tends to a constant value which 
does not depend on U but only on n according to the law 



respect to h and letting h going to zero: 



lim Xc = (^{n), 

1 ^OO 



(35) 



where a{n) — n{2 — n)/2. 



2. Spin susceptibility 

The spin magnetic susceptibility Xs can be computed 
by introducing an external magnetic field h, taking the 
derivative of the magnetization m = (ri|(j) — "•j,(i)) with 



Xs = 



dm 
'dh 



h=0 



The addition of a homogeneous magnetic field does not 
dramatically modify the framework of calculation given 
in Sec. |TT1 once one has taken into account the break- 
down of the spin rotational invariance. Details of the 
calculations will be presented elsewhere [61j: here we 
shall present only some results for the spin susceptibil- 
ity Xs relevant to our discussion. In Figs. I26b -c we plot 
the spin susceptibility as a function of the temperature 
for three representative values of C/ (C/ = —V, U = V 
and U = 2.1V^, respectively) and for different values of 
the filling (n = 0.25, 0.5, 0.75, 1). In the first region 
[U < 0, phase (a)], the spin susceptibihty vanishes at 
zero temperature for all values of the filling: at T = 
all electrons are paired and no alignment of the spin is 
possible. By increasing T, the thermal excitations break 
some of the doublons and a small magnetic field may in- 
duce a finite magnetization: Xs augments by increasing 
T up to a maximum, then decreases. Calculations show 
that, in the limit of vanishing temperature, the behavior 
of Xs in the phase (a) is 



lim Xs 



-0{f^-U) 



T 



The exponential behavior exhibited by the spin suscep- 
tibility is easily understood as a finite amount of energy 
is required to break the doublons. In the second region 
[0 < U < 2V, phases (&) and (c)], Xs diverges at T = 
for all values of the filling, but n = 1. At zero tempera- 
ture, for n < 1 all the electrons singly occupy a macro- 
scopic number of sites and their spins will align when 
the magnetic field is turned on. In the limit of vanishing 



18 



d" 




T=0. 1 


0.30 


V=l 

- U=-l _.-■-■-.. _^ 


T=0.5 

T=1.0 








0.20 
0.15 
0.10 
0,05 


-// '■■■■, 'v^ 





T=U.2 


T-l.O 


ij:::: 


^ 



(a) 



(b) 




(c) 



FIG. 27: The spin susceptibility as a function of the particle 
density for V = 1 and different values of the temperature, 
and for (a) U = ~1 ; (h) U = 1; (c) U = 2.1. 



temperature, 
filling, as 



lim y, = 



Xs diverges, according to the value of the 



for < n < 0.5 
for 0.5 < n < 1 



phase (&) 
phase (c). 



By increasing T, the alignment is disturbed by the ther- 
mal excitations and, as a consequence, Xs decreases. On 
the other hand, for n = 1 at T = 0, all electrons are 
paired and Xs = 0. At finite temperatures, Xs increases 
linearly up to a maximum, then decreases. In the third 
region [2V < U, phases (&) and (d)] at T = and for 
all values of the filling, all the sites are singly occupied: 
the spin susceptibility decreases with T and diverges at 
T = 0. It can be shown that in this region the spin 
susceptibility diverges as T — > according to the law 
limT^o Xs =^ n/T. It is easy to check that, for high tem- 
peratures, the spin susceptibility decreases with the Curie 
law Xs ~ T^^ in the entire (t/, n) plane: 



lim Xs 



a{n) 



where a{n) is the same [/-independent function appear- 
ing in Eq. (I35p . As a consequence, in the limit of high 
temperatures, the ratio Xc/Xs is a universal function of 
T: 



lim 



T. 



In Figs. [27k-c we plot the spin susceptibility as a function 
of the particle density, again for the three representative 
values oiU {U = ~V, U ^ V and U ^ 2. IF, respec- 
tively) and for different temperatures. In the first two 
regions Xs has a similar behavior: it increases with n up 
to a maximum around n = 0.5, then decreases, exhibit- 
ing a minimum at n = 1 , where at T = all electrons are 





FIG. 28: The density of states N{lu) as a function of the 
frequency uj at quarter filling, V = l, T = 0, 0.5, 1 and 2 for 
(a) l7=-l/2 ; (b) [7 = 1/2. 



paired. In the third region Xs is an increasing function 
of n and presents a maximum at n = 1 where, at T = 0, 
the maximum number of singly occupied sites is present. 



G. Density of states 

The behavior of the density of states is, obviously, 
richer at finite temperatures since the thermal excitations 
favor the emergence of new peaks. As one can clearly see 
from Fig. I28[ by increasing T the number of peaks in- 
creases leading to a larger number of states which the 
electrons can occupy. In fact, for T 7^ 0, the weights E 
appearing in Eq. (j28[) are all finite. Here we report, as 
an example, only the DOS at quarter filling in the neigh- 
borhood of [/ = 0. One obtains similar diagrams also for 
other values of the filling. 



H. Specific heat 



The specific heat is given by 

dE 



(36) 



where the internal energy E can be computed as the ther- 
mal average of the Hamiltonian ([T]) and it is given by 



E = UD + 2V\^^K 



(37) 



The specific heat exhibits a very rich structure in corre- 
spondence of the zero-temperature phase diagram shown 
in Fig. [TJ The possible excitations of the ground state 
are creation and annihilation of particles in either a singly 
or a doubly occupied site, induced by the Hubbard op- 
erators i}} and ■(/;'■''•', respectively. The corresponding 
transition energies are given in Eq. ([7]) and the peaks 
exhibited by the specific heat correspond to these tran- 
sitions. One may distinguish between the transitions in- 
duced by V'^^^ or by '0^''^ by looking at the position of the 
peaks, i.e., if the position changes or remains constant by 
varying U . The behavior of the specific heat at n = 0.25 
as a function of the temperature is shown in Figs. I29b -d. 
The specific heat presents one pronounced maximum for 
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FIG. 29: The specific iieat C as a function of tlie temperature 
T for y = 1, n = 0.25 and (a): U = -3,..., -0.5, (b): U = 
-0.5,..., -0.1; (c): U = 0.1,..., 0.5; (d): C/ = 1, ..., 5. 
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FIG. 30: The temperatures Ti and T2 at which the specific 
heat exhibits a peak are plotted as functions of U for n = 0.25. 



large negative values of L'^ at a temperature Ti (see Fig. 
E5k). The position of the maximum decreases linearly 
with U until [/ = 0. Around U ~ -~0.3V a second peak 
appears at a higher temperature (see Fig. [^5b). In this 
region, both peaks represent transitions mainly induced 
by ip^^\ The height of the first peak decreases and sud- 
denly goes to zero at ?7 = 0, where it vanishes. At U ~ 
there is only one peak located at T2 ~ 0, idV. For U > 
there is a small region where the double peak structure 
is again present (see Fig. [^5b). When U > 0.2F the first 
peak disappears while the position of the second peak 
slightly increases with f/ up to [/ « 2V where it assumes 
a constant value (see Fig. [53Ji) . The behavior of Ti and 
T2 as a function of U is shown in Fig. [3D1 The behavior 
of the specific heat for this value of the filling is easily 
understood if one recalls that at T = 0, by varying U, 
a FT is observed at [/ = from a state with doubly 
occupied sites to a state with singly occupied ones (see 
Figs. [3] and [8]). Accordingly, at finite temperatures, for 
U < the excitations are mainly induced by the Hub- 
bard operator -0^'') (Ti varies linearly), while for U > 
the excitations are mainly induced by ip {T2 is constant 
for U >2V). 

In Figs. [HTb-d the temperature dependence of the spe- 
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FIG. 31: The specific heat C as a function of the temperature 
T for V = 1, n = 0.5 and (a): U = -3, ..., -0.5; (b): 



(7 = -0.5,, 



-0.1; (c): [7 = 0.1,. ..,0.5; (d): 17 = 1,..,5. 



cific heat is shown for n = 0.5. The first peak Ti, appear- 
ing for large negative values of U (see Fig. 15Tb). decreases 
following a linear law up to C/ « —2V. At C/ ~ — 1.2T^, a 
second peak appears at a temperature T2 lower than T\. 
By decreasing |U|, T^ decreases linearly (see Fig. [5Tb V 
At C/ = there is only one peak situated at T\ « 0.56^. 
For small positive values of t/, the specific heat has again 
two peaks (see Fig. 15Tb). By increasing [/, T^ increases 
following a linear law up to [7 « 0.5V^; further increasing 
[/, T2 deviates from the linear behavior and tends to the 
constant value T2 = 0.21 for U > 2V (see Fig. |3lji). 
The position Ti of the first peak decreases by increasing 
U and, for U > 0.5^, the peak is no longer observable. 
The behavior of Ti and T2 as a function of U is shown 
in Fig. [32l It is worth noticing that, for positive U, the 
double peak structure is associated with the presence of 
a crossing point in the specific heat curves. A double 
peak structure, as well as the appearance of quasiuniver- 
sal crossing points, was already noticed in Refs. [l3, l63 | 
and in Sec. 3.4.2 of Ref. [5^ in the Hubbard model, in 
Ref. [i^ in an extended Hubbard-like model, and also in 
some experimental data j63||. When U > V/2 the dou- 
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FIG. 32: The temperatures Ti and T2 at which the specific 
heat exhibits a peak are plotted as functions of U for n = 0.5. 
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FIG. 33: The specific iieat C as a function of tlie temperature 
T for V = 1, n = 0.75 and (a): U = -3,. ..,-0.5; (b): U = 
-0.5,. ..,-0.1; (c): U = 0.1,. ..,0.5; (d): 17 = 1, 1.5, 1.75, 1.8 
and 2; (e): 17 = 5,. ..,10. 



ble peak structure and the crossing point disappear. As 
for the case n = 0.25, the elementary excitations which 
mainly contribute to C for n = 0.5 are those induced by 
V'^"^ for L/ < and by V-^*^ for C/ > 0. 

For n = 0.75, as it is shown in Figs. [55b-e. one ob- 
serves a double peak scenario both around t/ = and 
U — 21/, in correspondence of the two different transi- 
tions observed at T = 0. For negative values of C/, the 
specific heat has a peak at a temperature Ti (see Fig. 
155k ) which decreases following a linear law until [/ « 0. 
Around U k, —V, at a lower temperature T2 one finds a 
second peak which decreases following a linear law (see 
Fig. [33b V At U — there is only one peak situated at 
Ti w 0.62F. For small positive values of U the specific 
heat has again two peaks (see Fig. [SBbV By increasing 
U, T2 increases following a linear law up to C/ « 0.35^, 
where this peak is no longer observable. On the other 
hand, the position of the first peak Ti decreases by in- 
creasing U wp to U = 2V , where this peak disappears. 
When U approaches the value U — 2V a third peak ap- 
pears at some temperature T3 (see Fig. [33H). We recall 
that at C/ = 2y and r = (see Figs. M and ^ there 
is a transition from an ordered state characterized by al- 
ternating empty and occupied sites to a disordered state 
with only singly occupied sites. T3 starts at t/ ~ 1.75^ 
and increases linearly with U. Further increasing the 
on-site potential one finds that at f/ w 5AV a fourth 
peak, situated at a lower temperature T4, appears (see 
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FIG. 34: The temperatures Ti, T2, T3 and T4 at which the 
specific heat exhibits a peak are plotted as functions of U for 
n = 0.75. 



Fig. [33b). As a result, in the region U > 5AV the spe- 
cific heat exhibits a double peak structure, with T3 in- 
creasing linearly while T4 remains constant at the value 
T4 « 0.4iy. In Fig. [Mlwe show the behavior of Ti, T2, 
T3, T4 as functions oi U . In contrast to what happens in 
the range < n < 0.5, for n — 0.75 the excitations in- 
duced by ■)/;*•'''' and by ■i/'^^-' are both present for U > 2V: 
indeed, at T = 0, one observes both singly and doubly 
occupied sites. The behavior of the specific heat at n = 1 
as a function of temperature is shown in Figs. 135b . and 
135b. When T = and [/ = 2V, there is a FT from an 
ordered state with doubly occupied sites organized on a 
checkerboard distribution to a homogeneous state, where 
all the sites are singly occupied. The possible excitations 
are creation and annihilation of doubly occupied states. 
One observes a peak for negative values oi U at a. tem- 
perature Tl, as shown in Fig. 136b . Ti decreases linearly 
until [/ w 0. As U gets close to 2V, Ti further decreases, 
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FIG. 35: (a) The specific heat C as a function of the temper- 
ature r for V = 1, n = 1 and (a) U = -2, 0, 1, 2, and 4; (b) 
U = 1.8, 1.82, 1.84 and 1.86. 

while its height hi exponentially increases and vanishes 
at [/ = 2V. T2 appears at C/ « 1.81^ and increases 
linearly with U. When U > 2V, further increasing the 
on-site potential leads to an increase of the thermal en- 
ergy required to excite an rj electron above the ground 
state. Thus, the peak in the specific heat shifts towards 
higher temperatures as it is clear from Figs. 135b and l36b . 
In Figs. [36k-b we report the behavior of Ti, T2 and hi 
as a function of U. Both Ti and T2 vary linearly with U, 
according to the fact that the only possible excitations 
are those induced by -0^''^. 

The divergence of the specific heat in the limit U -^ 
2V, when T — > 0"*", hints at a first order phase tran- 
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(a) (b) 

FIG. 36: (a) The temperature Ti and T2 as functions of the 
on-site potential for T^ = 1, n = 1. (b) The height hi as a 
function of the on-site potential for V — 1, n = 1. 



sition. Another way of looking at this is to study the 
behavior of the internal energy. In Fig. [37] we plot the 
internal energy E as a function of T for several values 
of U just below U — 2V. One sees that E exhibits a 
jump around the related Ti temperatures. By increas- 
ing U, the jump becomes steeper and, at the same time, 
smaller. In the limit U — > 2V, E exhibits a discontinuity 
with Ti — > and the specific heat diverges. On the other 
hand, at U = 2V, the specific heat is finite with a single 
peak. Sharp peaks in the specific heat at half filling, in 
the neighborhood of the phase transition, have also been 
observed in Refs. [2J, I44I, although with approximate 
methods. 

To end this subsection, we want now to comment on 
the issue related to the two-peak structure exhibited by 
the specific heat. This issue has been largely discussed 
in the literature in the context of the ordinary Hubbard 
model, where the common interpretation explains the 
two peaks as due to spin and charge excitations. This 
explanation does not apply to the AL-EHM where only 
charge excitations are present. Our exact solution of the 
model points out that the possible excitations [recall the 
energy spectrum given in Eq. ([7])] are due to a redistri- 
bution of the charge density induced by the creation and 
annihilation of electrons under the constraint imposed 
by the preexisting charge distribution. For instance, the 
operator ^^(i)n"(z) creates an electron at the site i, pro- 
vided that this site is empty and the neighbor site(s) is 
(are) occupied. 

Our extensive study of the U dependence of the tem- 
perature maxima T„ (at which the specific heat exhibits 
a peak) clearly indicates that the observed peak struc- 
ture is due to charge excitations induced by the two op- 
erators ^'■^^ and ■0'^''^. An explanation of the two peak 
structure, as due solel y to charge excitations, has also 
been proposed in Ref. [l^ , where an exact solution of a 
Hubbard-like model has been given. 



I. The entropy 

The standard way to compute the entropy is via the 
derivative of the specific heat: 
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FIG. 37: The internal energy E as a. function of T for V — 1, 
n — 1 and various values of U. 



This expression requires the knowledge of the zero tem- 
perature entropy S(0), which can be computed by means 
of the formula S{0) = feeing, where ks is the Boltz- 
mann's constant and g is the number of accessible mi- 
croscopic states, i.e., all the states with energy equal to 
the ground state energy. Now, if the lowest energy level 
is non degenerate, g = I and 5(0) = 0. On the other 
hand, if the lowest energy level is degenerate g ^ I and 
S{0) y^ 0. This latter situation is the one occurring in 
the AL-EHM: as it is shown in Figs. H E and [12l at 
T = 0, every phase (except for n = 1 and U < 2V) ex- 
hibits a macroscopic degeneracy. As a result, the zero 
temperature entropy is different from zero. One state- 
ment of the third law of thermodynamics is that, as the 
temperature approaches zero, the entropy becomes inde- 
pendent of the external parameters. This still holds for 
the AL-EHM confined in one of the possible phases. As 
we will show later, S(0) depends on the particle density 
and on the on-site potential since, by varying n or U , the 
system can undergo a PT. 

In some special cases, it easy to compute S{Q). For 
example, from Figs. [3l [8l and [T2l it is easy to see that: 

(i) for U > and n = 0.5, the ground state corresponds 
to a checkerboard distribution of singly occupied sites 
with arbitrary spin; the ground state energy is indepen- 
dent on the spin orientation and the number of accessible 
states is 5 = (2)^/2; therefore S{0) = f ln(2). 

(ii) for U < 2V and n — 1, the ground state corre- 
sponds to a checkerboard distribution of doubly occupied 
sites. There is no degeneracy and, as a result, the number 
of accessible states is 5 = 1 and 5(0) = 0. 

(iii) for U > 2V and n — 1, the ground state corre- 
sponds to an homogeneous distribution of singly occu- 
pied sites with arbitrary spin; the ground state energy is 
independent on the spin orientation and the number of 
accessible states is g — (2)^; therefore 5(0) = 7Vln(2). 

Discarding few simple cases, it is generally not easy to 
compute 5(0). As a consequence, one has to use different 
formulas to compute the entropy. In the limit of infinite 
temperature, the entropy can be computed as S{oo) = 
ks In goo, where goo is the total number of states. For the 
model under consideration it is easy to show that 

lim 5(T) = 21og2-nlog(n) -(2- n)log(2-n). (39) 
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FIG. 38: The entropy 5* as a function of the temperature T 
for V = 1, [/ = -1, 1, 2.1 and (a): n = 0.25, (b): n = 0.5, 
(c): n = 0.75, (d): n = 1. 



Here and in the following we set the Boltzmann's con- 
stant ks = 1 and consider the entropy per site. For 
general values of the temperature, recalling that (n, — /i) 
are thermodynamically conju gat ed variables, one can ex- 
press the entropy per site as [59| 



Sin,T,U)^- f 
Jo 



7:7;:. dn 

dT 



1 
T 



E{n,T,U)- I fi{n',T,U)dn' 




(40) 



This expression requires only the knowledge of the in- 
ternal energy and of the chemical potential. By making 
use of Eq. (l40l) . we shall present in the following sev- 
eral results for the entropy as a function of T, n, and 
U. In Figs. [38] we plot the entropy as a function of the 
temperature for relevant values of the particle density n. 
One clearly sees that S is always a decreasing function 
of the temperature and that, for n < I, S does not van- 
ish for T ^ 0. For instance, for n = 0.75 at T = 
three phases are present by varying U: accordingly, the 
entropy for T — > tends to three different constants (see 
Fig. [55b ) which depend only on the degeneracy of the 
related ground state. Furthermore, it is easy to check 
that, in the limit T — > 00, the entropy tends to the value 
given in Eq. ([5^ . It is interesting to note, as it is shown 
in Fig. [351 that at half filling, tuning the parameter U in 
the neighborhood of the critical value U = 2V, the en- 
tropy suddenly vanishes at the same temperature Ti at 
which the specific heat increases exponentially (see Fig. 
[55b). This can be easily understood by recalling that 
C — T ■ dS/dT: the divergence of the specific heat is a 
consequence of the discontinuity of the entropy. One can 
think of this divergence as induced by a change from a 
non-degenerate state to an infinitely degenerate state. 
When plotted as a function of the particle density (see 




FIG. 39: The entropy as a function of the temperature for 
V — 1, n — 1 and different values of U close to 2V . 



Fig. I40p . the entropy shows a maximum around quarter 
filling when one considers an attractive on-site interac- 
tion. At low temperatures, the entropy increases by in- 
creasing the particle density; it has a maximum around 
quarter filling and then decreases vanishing, as expected, 
at half filling. On the other hand, considering a repul- 
sive on-site potential, the entropy presents a minimum 
at quarter filling and two maxima around n — 0.25 and 
n = 0.75. As evidenced in Fig. [iOb . when U > 2V the 
entropy does not vanish at half filling. The increasing and 
decreasing of the entropy by varying the particle density 
can be understood as due to the increase and decrease of 
the number of accessible states, as it is easy to check by 
looking at Figs. [3l [H and [T2| It is worthwhile to observe 
that the minimum values of the entropy are observed in 
correspondence of the ordered states, where the number 
of accessible states is, of course, minimum. 





(a) 



FIG. 40: The entropy S as a function of the particle density 
n for y = 1 and r = 0.1 (a): U = -3, -1, -1/2, and 0; (b): 
U = 0, 1/2, 1, 3/2 and 5/2. 

A useful representation of the finite-temperature phase 
diagram is obtained by plotting the entropy as a func- 
tion of the on-site potential U, as it is shown in Figs. 
[4Th-d. The U dependence of the entropy is rather dra- 
matic in the neighborhood of the values at which a zero 
temperature transition occurs. At low temperatures, the 
entropy presents a steplike behavior and becomes rather 
insensitive to variations in U for sufhciently large on-site 
repulsive and attractive interactions. 

When n < 0.5, the entropy presents a peak around 
U — which becomes more pronounced as the tempera- 
ture decreases. At quarter filling, one observes a decrease 
of the entropy by varying U since, in the limit T — > 0, the 
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FIG. 41: The entropy 5 as a function of U ior V — 1, for 
temperatures varying in the interval (0,2) and (a): n — 0.25, 
(b): n = 0.5, (c): n = 0.75, (d): n = 1. 



electrons are organized in an ordered checkerboard distri- 
bution of singly occupied sites. In the region 0.5 < n < 1, 
at low temperatures, one finds three regions of constant 
entropy corresponding to the three phases present at 
T = when one moves, at fixed particle density, from 
negative to positive values of U. At half filling, as ex- 
pected, one finds a rather sharp increase of the entropy 
only near U = 2V. In Figs. [Uwe report the value S{0) 
obtained by subtracting Eqs. p8|) and (|40l) . The study 
of 5(0) constitutes another way to detect the zero tem- 
perature transition from thermodynamic data. In fact, 
as it is clearly shown in Fig. 311 for each value of the 
electronic filling, the zero temperature entropy presents 
a discontinuity right at the values of the on-site potential 
at which one observes a phase transition. 



V. CONCLUDING REMARKS 

We have evidenced how the use of the Green's function 
and equation of motion formalism leads to the exact solu- 
tion of the one-dimensional extended Hubbard model in 
the atomic limit. We provided a comprehensive and sys- 
tematic analysis of the model by considering all the rele- 
vant correlation functions and thermodynamic quantities 
in the whole space of parameters n, T, and U (having 
chosen F = 1 as the unity of energy), including attrac- 
tive on-site interactions. This study has shown that, at 
zero temperature, the model exhibits phase transitions 
for specific values of n and U . In particular, we have 
identified four phases in the (t/, n) plane and PTs are ob- 
served at the borders of these phases. Various types of 
long-range charge ordered states have been observed: (i) 
at half filling and U <2V , a checkerboard distribution of 
doubly occupied sites; (ii) at quarter filling and C/ > 0, 
a checkerboard distribution of arbitrary-spin singly oc- 
cupied sites; (iii) in the entire region 0.5 < n < 1 and 



< t/ < 2y, an ordered state with alternating empty 
and occupied sites. The onset of a CO is signalled by the 
divergence of the correlation length in the charge corre- 
lation function. It is worthwhile to mention that the PT 
observed at half filling and U -^ 2V \s a, first order tran- 
sition, as confirmed by the divergence of the specific heat 
and by the discontinuity in the entropy. 

We derived the phase diagram by analyzing at T = 
various local correlators which provide information on 
how the electrons are distributed on the sites of the chain. 
As shown in details in Sees. Ill and IV, the phase di- 
agram and the observed CO states can be equivalently 
studied by analyzing the behavior at T = of other quan- 
tities. The chemical potential and the double occupancy 
show discontinuities where PTs occur. The existence of 
a long-range order is well evidenced by the periodicity of 
the charge and double occupancy correlation functions. 
Another way to detect the zero temperature phase dia- 
gram is provided by the study of the behavior in the limit 
T — ^ of the charge and spin susceptibilities and of the 
entropy. In particular, the entropy at zero temperature 
is non zero (except for n = 1 and U < 2V) and is deter- 
mined only by the degeneracy of the ground state. In the 
limit T ^ the spin susceptibility diverges exponentially 
or with a power law, depending on the phase. 

In Sec. IV we have studied the model at non zero tem- 
peratures. At finite but low temperatures the properties 
of the system are strongly determined by the zero tem- 
perature solution: a finite range order persists for a wide 
range of T, as evidenced by the behavior of the correla- 
tion functions. As a function of n, the charge suscepti- 
bility and the entropy exhibit a two-peak structure with 
minima at n = 0.5 and n ~ 1, where PTs are observed. 

The use of equation of motion and the introduction of 
a closed set of composite operators ip^^^ and ip^'^\ eigen- 
operators oi H, allowed for an exact determination of the 

and of the relative 
As a result, the density of states 
has been exactly computed for any temperature, and we 
performed a detailed analysis of the specific heat C In 
particular, by studying the specific heat in the parameter 
regions corresponding to the four T = phases, it has 
been possible to identify the excitations contributing to 
C as due to a redistribution of the charge density induced 
by the operators "0 and V' through the creation and 
annihilation of electrons, provided that a charge distri- 
bution preexisted. The presence of a two-peak structure, 
observed in some regions oi U, is explained as due to 
charge excitations induced by different components of the 
multiplet field operator -0(''). 
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APPENDIX A: THE PARAMETERS k^*^^ AND A^''^ 

As already mentioned in Sec. [TTl the use of the formal- 
ism of Green's functions and equations of motion requires 
two ingredients: composite operators and algebraic prop- 
erties. In this Appendix we shall summarize the main 
results which can be obtained by exploiting the algebra 
satisfied by the operators. By means of the algebraic 
rules 

[Dii)]" ^ Dit), (Al) 

where ak — 1^ — 2 for A; > 1, it is possible to derive the 
recursion rule 



'W]'-E^™[^"«r 



(A2) 



Ak) 



where the coefficients Am are rational numbers, given 
by 



iW 



A^^'^ = ^ i,J -104 + 57 ■ 2'=+! - 56 ■ 3*^ + 11 • 4'=) 



3 • 2^+1 



A 



(fc) 



1 



i(fc) 



3-2^ 
1 



— (18 - 3 • 2''^+^ + 14 ■ 3*^ - 3 • 4'=) 



3-2' 



— (-4 + 3-2''+^ -4-3'= + 4*^) 



(A3) 



It is easy to show that the relation X]m=i ^™ ~ ^ holds. 
As it has been shown in Sec. IlICi by using the i/o- 
representation, the expectation value of a generic opera- 
tor O can be expressed as 



(0) = 






-mi 



(A4) 



Recalling that if/ = 2Vn(%)n°'(i) and by using the recur- 
sion rule (IA2I) one obtains 



(A5) 



m— 1 



where 






oo 



Y^{-iy^la^A(-\2pvr 



(A6) 



By recalling Eq. (|A3p . one immediately sees that the 
coefficients /,„ and gm are polynomials of i^ = e"'^^ and 



n=2 



they are explicitly given by 
1 
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, 28 . o 10 

/3 = -2K'' + —K^ - 16K^ + UK - —, 

2 A o o o o 2 

^^-3^ -3^ +^^-3^+3' 



(A7) 



and 



31 = — (K - 1)-V25 + 29K + 15K^ + 3K^), 
6 

c,2 = l(-35 + 208ii: - 332K^ + 112K^ 
6 

4 I Qi^5\ 

4(' 1 I ,1 J^ I « 1^2 I ^ T^-S I 7^4\ 



-^927^"*- 56X^-1- 11 is:*), 

2 
33 = --{K ~ 1)3(5 - 21K - \2K'^ + 4/s: 

+ 9i<i:'* + 3ii'^), 

2 



(A8) 



g^= [K - 1)'*(-1 -I- 4ivr + 6is:^ + 4^3 + K 

o 

One can now easily compute the parameters k,^^' and 
A^*^). By using Eqs. (|A4p and (jA5|) and recalling the 
properties of the Hq representation |54| one has 



p=l 



A('^)-:^{Bi(K(*)f)o 

4 

+ E[(5i + 2i?2)/p + 2i32.gJ(K(i)]''+'=)o}, 



(A9) 



where 

4 

To = (e-^^Oo = 1 + E(^i/p + 525p)K(*)]^ 

„ , ,.,, 26^5^(1 -He/'f^-^)) (AlO) 



^2 = (i?(*))0 



o/3(2m-C/) 



1 + 2e/3A' + e/3(2A'-tf) ■ 



By using the properties of the Hq representation, the 
expectation values ([n"(z)]'^)o can be written as 

(K(^)]3)o = \x, + \x^ + \xl + \x,X^, 

(K(^r)o = ^^i(l + IX,) + Ja2 + \x,X2 + ^x|, 

(All) 
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where the two parameters Xi and X2 are defined as 

(A12) 



Xi = (n"(z))o, 
X2 = {D''i^))o■ 



purposes, correlation functions (|B1|) can be written as: 

3 



Then, one has 



To = l-Bi+B2 + iBi- 2B2){1 + aXi + a^Xa)^ 
+ B2{l + dXi+d^X2f, 

with a = K-1 d^K'^-l. The correlators {[n°'{i)f)Q 
with fc > 4 can be computed by using the recursion rule 
(IA2[) . The previous formulas show that the local CFs 
k'^^^ and A''^-' are exactly known in terms of the two basic 
parameters Xi and X2. A similar procedure can be used 
for computing the other parameters tt^*^^ = ([Z?"(i)]'^), 
JC^) = {n{i)[D''{i)f)/2, 0'-''^ = {D{i)[D°' {{)]'') /2. In par- 
ticular, one obtains 



B2X2 
2To 



(1 + 4a + 46 + c)[l + (2a + &)Xi + (c + 26)^2 



^(2) = 1^(4) _ ^(3) + ii^(2) _ i„ ^ i^ 

3 12 4 2' 

^ -{Si(l + 2a + fe) 



^^^^ = ^M.)[z^"W]) = 2^^_,,,^^ 



• (1 + aXi + bX2) + 2B2[6a^Xi + fe^(4Xi + 7X2) 
+ c(l + X2 + CX2) + a(2 + Xi + II6X1 + 2cXi 
+ 66X2 + 4cX2) + 6(3 + Xi+ cXi +X2 + 6CX2)]}, 

(A13) 

where b = {K ~ l)^ and c = 4is: - 4K^ + K'^ - 1. 



APPENDIX B: NONLOCAL CORRELATION 
FUNCTIONS 

One can define the following nonlocal correlation func- 
tions: 



e(o 

p(l 



*,i) = {n{i)n{j)) == ^a^y£, 
3 

1=1 

3 
*,j)-K(*)ri(j))=^d,y,, 
f=l 

3 

e=i 
3 
*,j)-(i?"(*)"(j))=^m,r,, 
i=l 

3 

*,j) = (i?(»)i?(j)) = E^^^^' 

£=1 
3 

z,j) = (i^(j))=^QZ,, 
3 

3 

1=1 
3 



(B2) 



A«(j,j> 


= (n(*)K(i)]'=n(j)), 


K^^\i,i)- 


= {[n^{z)fn{3)), 


e(^)(*,j) = 


- {D{i)[D'^{i)fn{j)), 


nW(z,j). 


- {[D-{^}]^nm, 


mW(z,j). 


- {D{i)[D'^{i)fD{j)), 


pW(z,j) = 


- {[D'^{t)fD{j)), 


Q('H^,j) = 


- {n{t)K{t)fD{j)), 


i?(^)(z,j). 


- {[n-{^tD{J)), 



where, by assuming unitary lattice spacing, the function 
Yi and Zi are given by 

Yi = {n{j))o, Zi - {D{j))o, 

Y2 = (n(^ + l)n(j))o, Z2 - {n{i + l)D{j))o, (B3) 

Y3 = {Dii + l)n(i))o, Z3 = {D{i + l)DiJ))o. 

The coefficients a^, 6^, q, d^ and nii appearing in Eqs. 
(|B2p can be explicitly computed by means of the formu- 
las (|A4p . (jA10[) and by using the properties of the Hq- 
representation. As an example, the coefficients a£ are 
given by: 



ai 



Bi{l + aXi + bX2) + 2B2 [(a + b)Xi + {b + 0)^2] 



(Bl) a2 



as 



where j is a site m-steps (to > 1) distant from the site 
i. By using the iJo-representation, the relevant, for our 



To 
a2(Bi + 6B2)Ai + 26B2[1 + CA2 + 6(Ai + 2A2)] 


To 
a[Bi(l + 6A2) + 2^2(1 + 46Ai 


+ 36A2 + 2CA2)] 


To 
62(4B2Ai + S1A2 + 6S2A2) + 2cB2{l + 2aXi + 0X2) 


To 
b[Bi{l + aAi) + 2^2(1 + 3aAi 


+ cXi + 4cA2)] 



To 



(B4) 
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where the coefScients a, b and c have been defined in 
Appendix [XI To determine the functions Yi, Y2, Y^ and 
Zi, Z2, Zj, one requires translational invariance which 
leads to the foUowing equations: 



(r 


iW) = 


' Hj)) 


=^ 


n = 


>> 

£=1 


lYu 


(^ 


«) = 


{D{j)) 


=> 


D = 


3 

= > ]c,Z,, 
1=1 


(i^(z)n(j)) 


= W 


^)D{j)) 


=> 


3 

>; 
3 


biYi^ 


3 

-- ^ ' Of Zf , 
f=i 
3 


"(0^(j)) = 


= (n"( 


^)D{j)) 


=> 


>: 

1=1 


mtYt 


(B5) 



The charge and double occupancy correlation functions 
can more conveniently written in terms of the relative 
distance m = \i — j\ as 



A(°)(to) = {n{i)n{i + m)) = ^ mYi, 

1=1 
3 
M^^\m) ^ {D{i)D{i + m)) =Y^ b^Zi 



(B6) 



1=1 



By solving the system (|B5l) - (jB6l) one finds 

Yi ^uen + ViD + WiA^°\m) + zi,M^°\m), 
Zi^pin + qeD + r^ A(") (m) + s^M^o) (m) , 



(B7) 



where the coefficients M£, u^, w^, zi and p^ , g^ , r^, si are 
just functions of the coefficients Of, 6^, q, di and to^. 
Next, we note that by definition 



K^^^i,j) ^- [A(")(m + 1) + A(")(TO-1) 
P(i) (i, j) = - \m'-°^ (m + 1) + M(") (m - 1) 



(B8) 



Upon substituting Eqs. ([62]) and ([BT]) in Eq. (|B8]) one 

gets 

A(°)(m + 1)+A(°)(m-1) 

3 
= 2^di\uen + ViD + weA^°^ (m) + z^A-f (°) (m) 

Af(°)(TO + l) + A-/(°)(m-l) 
3 



e=i 



(B9) 



These equations give the sought recurrence relations for 
the correlation function A(°)(m) and M(*'^(r7i). As a re- 



sult, one can write these equations as 

A(") (m + 1) = VKo A(") (to) - A(") (to - 1) 

M(") (to + 1) = OoAf (") (to) - Af (") (to - 1) 
+ eiA<^°\m) + 02, 

where one has defined 

3 3 

Wq = 2'^dewe, eo = 2^TOfS^, 



(BIO) 



Wi=2Y,dtze, 



£=1 
3 



61 = 2^mgri, 

f=i 
3 



W2 ^ 2^^de[uin + vi D], 82 = 2^^mi[pin + qi D]. 
e=i 1=1 

(Bll) 

Now, upon defining 

G(m)=A(°)(TO)-n2, 



F(m) ^ M^°\m) - D^ , 

one can rewrite Eq. (jBlOp as 

G{m + 1) = WoG{m) ~ G{m - 1) + WiF{m) 
+ W2 + Wqt? - 2r? + l^iL'^ 

i^(TO + 1) = eoi^(TO) - i^(TO - 1) + eiG(m) 



(B12) 



(B13) 



It is not difficult to show that the following identities do 
hold for any n: 



W2 + Wan^ - 2r? + WiD'^ = 0, 
62 + 60!?^ - 2L»2 + Oin^ = 0. 

Therefore, Eqs. (JB13P simplify as 

G(to + 1) = WoG(m) - G{m - 1) + WiF{m), 
F{m + 1) = eoF(m) - F{m - 1) + eiG(m). 



(B14) 



(B15) 



The results (jB14p are just a manifestation of the ergodic 
conditions 

lim A^ -^(to.) = lim {n{i)n{i + rn)) ~ {n{i)) , 

m — >oo m — >oo 

lim Af(")(TO)= lim {D{i)D{i + m)) ^ {D{i)f , 

rn —>-OG 7n—>-OG 

(BI6) 

which break down when there is a long range order. Ex- 
amples of this scenario are shown in Sec. IIIII Equations 
(JBISP together with Eq. pT|) allow one to compute the 
correlation functions A*^°)(m) and M^^^{m) for any value 
of TO,. Then, by means of Eq. (|B7p one can compute the 
functions Yi, Y2, I3, and Zi, Z2, Z3, whose knowledge 
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allows one to compute any correlation functions, by us- 
ing Eq. (|B2p or relative formulas. Equations (|B15p show 
that G{m) and F{m) have the following form 



F{Tn):^Cp"' + Eq"', 



(B17) 



where the coefficients A, B, C, E, p, and q are given by 
Gq{F2Gq — F0G2) + 2Gi{FoGi — FiGo) + GoQ 



A = 



2Q 



B = Go-A, 

Fo{F2Go — F0G2) + 2Fi{F()Gi — FiGa) + FqQ 



G = 

E = Fo-C, 

F2G0 — F0G2 — Q 



2Q 



P = 



2(FiGo - FoGi) ' 

F2G0 — FqG2 + Q 

2(FiGo - FoGi) ' 



with 



and 



Q = [iF2Go - FoG2f 

+ 4(FiGo - FoGi)(FiG2 - F2G1)] '^\ 

Go = G(0) ^n + 2D-n^, 

Gi = G(1)==2A(i' -n^, 

G2 = G(2) = 2k(2) - n - 21) - n^, 



(B18) 



(B19) 



and 



Fi =F(1) = 26l(^) -D^^ 

F2 = F(2) = 27r(2) -D-D"^. 



One can rewrite Eq. (IB17[) as 



G(m) = A [sign(p)]" e 
F{m) = G [sign(p)]" e' 






i3[sign(g)]"e-'»/«°, 

■i?[sign(g)]"e-™/«-, 

(B20) 



where the two correlation lengths ^„ and S,d are defined 
as 



6. 



In 


-11 


"in 


-11 



(B21) 



As it has been shown in Sec. IIIIl in the limit of zero 
temperature the correlation length ^„ diverges in some 
regions of the (C/, n) plane, signalling the presence of a 
phase with long-range order. 
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